source:
issm/oecreview/Archive/25834-26739/ISSM-26664-26665.diff
Last change on this file was 26740, checked in by , 3 years ago | |
---|---|
File size: 15.7 KB |
-
../trunk-jpl/test/NightlyRun/test543.m
Cannot display: file marked as a binary type. svn:mime-type = application/octet-stream
1 %Test Name: PigTran Frontalforcingsrignotautoregression1 %Test Name: PigTranStochasticforcing 2 2 md=triangle(model(),'../Exp/Pig.exp',10000.); 3 3 md=setmask(md,'../Exp/PigShelves.exp','../Exp/PigIslands.exp'); 4 4 md=parameterize(md,'../Par/Pig.par'); … … 5 5 md=setflowequation(md,'SSA','all'); 6 6 md.timestepping.start_time = 0; 7 7 md.timestepping.time_step = 1; 8 md.timestepping.final_time = 25;8 md.timestepping.final_time = 10; 9 9 10 % separate domain in 2 basins11 idb asin= zeros(md.mesh.numberofelements,1);10 %Basin separation TF 11 idb_tf = zeros(md.mesh.numberofelements,1); 12 12 iid1 = find(md.mesh.x<=-1.6e6); 13 13 for ii=1:md.mesh.numberofelements 14 14 for vertex=1:3 15 15 if any(iid1==md.mesh.elements(ii,vertex)) %one vertex in basin 1 16 idb asin(ii) = 1;16 idb_tf(ii) = 1; 17 17 end 18 18 end 19 if idb asin(ii)==0 %no vertex was found in basin 120 idb asin(ii) = 2;19 if idb_tf(ii)==0 %no vertex was found in basin 1 20 idb_tf(ii) = 2; 21 21 end 22 22 end 23 % Basin separation default 24 idb_df = zeros(md.mesh.numberofelements,1); 25 iid1 = find(md.mesh.x<=-1.62e6); 26 for ii=1:md.mesh.numberofelements 27 for vertex=1:3 28 if any(iid1==md.mesh.elements(ii,vertex)) %one vertex in basin 1 29 idb_df(ii) = 1; 30 end 31 end 32 if idb_df(ii)==0 %no vertex was found in basin 1 33 idb_df(ii) = 2; 34 end 35 end 36 % Dimensionalities 37 nb_tf = 2; 38 nb_clv = 2; 39 nb_flmlt = 2; 23 40 24 % calving and frontalforcingsparameters41 %Calving parameters 25 42 md.mask.ice_levelset = 1e4*(md.mask.ice_levelset + 0.5); 26 md.calving=calvingvonmises(); 43 md.calving.calvingrate = 0.3*ones(md.mesh.numberofvertices,1); 44 md.levelset.spclevelset = NaN(md.mesh.numberofvertices,1); 45 md.levelset.migration_max = 10.0; %avoid fast advance/retreat of the front 46 %Frontal forcing parameters 27 47 md.frontalforcings=frontalforcingsrignotautoregression(); 28 md.frontalforcings.num_basins = 2;29 md.frontalforcings.basin_id = idbasin;48 md.frontalforcings.num_basins = nb_tf; 49 md.frontalforcings.basin_id = idb_tf; 30 50 md.frontalforcings.subglacial_discharge = 0.1*ones(md.mesh.numberofvertices,1); 31 51 md.frontalforcings.beta0 = [0.05,0.01]; %intercept values of TF in basins [C] 32 52 md.frontalforcings.beta1 = [0.001,0.0001]; %trend values of TF in basins [C/yr] … … 35 55 md.frontalforcings.ar_timestep = 2; %timestep of the autoregressive model [yr] 36 56 md.frontalforcings.phi = [[0.1,-0.1,0.01,-0.01];[0.2,-0.2,0.1,0.0]]; %autoregressive parameters 37 57 58 % Floating Ice Melt parameters 59 md.basalforcings.floatingice_melting_rate = 0.1*ones(md.mesh.numberofvertices,1); 60 61 % Covariance matrix 62 covtf = 1e-4*eye(nb_tf); 63 covclv = 1e-1*eye(nb_clv); 64 covclv(1,1) = 1/10*covclv(1,1); 65 covflmlt = 0.05*eye(nb_flmlt); 66 covglob = blkdiag(covtf,covclv,covflmlt); 67 38 68 %Stochastic forcing 39 69 md.stochasticforcing.isstochasticforcing = 1; 40 md.stochasticforcing.fields = [{'FrontalForcingsRignotAutoregression'}]; 41 md.stochasticforcing.covariance = 1e-4*[[1.5,0.5];[0.5,0.4]]; %global covariance among- and between-fields 70 md.stochasticforcing.fields = [{'FrontalForcingsRignotAutoregression'},{'DefaultCalving'},{'FloatingMeltRate'}]; 71 md.stochasticforcing.defaultdimension = 2; 72 md.stochasticforcing.default_id = idb_df; 73 md.stochasticforcing.covariance = covglob; %global covariance among- and between-fields 42 74 md.stochasticforcing.randomflag = 0; %determines true/false randomness 43 75 44 md.transient.ismovingfront = 1; 45 md.levelset.spclevelset = NaN(md.mesh.numberofvertices,1); 46 md.levelset.migration_max = 1e10; 76 md.transient.ismovingfront = 1; 77 md.transient.isgroundingline = 1; 47 78 48 md.transient.requested_outputs = {'default','Calving Meltingrate'};79 md.transient.requested_outputs = {'default','CalvingCalvingrate','CalvingMeltingrate','BasalforcingsFloatingiceMeltingRate'}; 49 80 md.cluster=generic('name',oshostname(),'np',2); 50 81 md=solve(md,'Transient'); 51 82 52 83 %Fields and tolerances to track changes 53 84 field_names ={... 54 'Vx1' ,'Vy1' ,'Vel1' ,' Pressure1' ,'Bed1' ,'Surface1' ,'Thickness1' ,'MaskIceLevelset1', 'CalvingMeltingrate1',...55 'Vx2' ,'Vy2' ,'Vel2' ,' Pressure2' ,'Bed2' ,'Surface2' ,'Thickness2' ,'MaskIceLevelset2', 'CalvingMeltingrate2',...56 'Vx10','Vy10','Vel10',' Pressure10','Bed10','Surface10','Thickness10','MaskIceLevelset10', 'CalvingMeltingrate10',...85 'Vx1' ,'Vy1' ,'Vel1' ,'Thickness1' ,'MaskIceLevelset1' ,'CalvingCalvingrate1' ,'CalvingMeltingrate1' ,'BasalforcingsFloatingiceMeltingRate1',... 86 'Vx2' ,'Vy2' ,'Vel2' ,'Thickness2' ,'MaskIceLevelset2' ,'CalvingCalvingrate1' ,'CalvingMeltingrate2' ,'BasalforcingsFloatingiceMeltingRate2',... 87 'Vx10','Vy10','Vel10','Thickness10','MaskIceLevelset10','CalvingCalvingrate10','CalvingMeltingrate10','BasalforcingsFloatingiceMeltingRate10',... 57 88 }; 58 89 field_tolerances={... 59 1e-11,2e-11,2e-11,1e-1 2,2e-11,6e-12,9e-12,1e-12,1e-9,...60 2e-11,1e-11,1e-11,9e-1 2,2e-1,2e-11,2e-11,1e-11,1e-9,...61 2e-10,1e-10,1e-10,1e-10, 1e-10,1e-10,1e-10,1e-9,1e-9,...90 1e-11,2e-11,2e-11,1e-11,1e-9,1e-10,1e-10,1e-10,... 91 2e-11,1e-11,1e-11,9e-11,2e-9,1e-10,1e-10,1e-10,... 92 2e-10,1e-10,1e-10,1e-10,5e-9,1e-10,1e-10,1e-10,... 62 93 }; 63 94 field_values={... 64 95 (md.results.TransientSolution(1).Vx),... 65 96 (md.results.TransientSolution(1).Vy),... 66 97 (md.results.TransientSolution(1).Vel),... 67 (md.results.TransientSolution(1).Pressure),...68 (md.results.TransientSolution(1).Base),...69 (md.results.TransientSolution(1).Surface),...70 98 (md.results.TransientSolution(1).Thickness),... 71 99 (md.results.TransientSolution(1).MaskIceLevelset),... 100 (md.results.TransientSolution(1).CalvingCalvingrate),... 72 101 (md.results.TransientSolution(1).CalvingMeltingrate),... 73 (md.results.TransientSolution(2).Vx),... 74 (md.results.TransientSolution(2).Vy),... 75 (md.results.TransientSolution(2).Vel),... 76 (md.results.TransientSolution(2).Pressure),... 77 (md.results.TransientSolution(2).Base),... 78 (md.results.TransientSolution(2).Surface),... 79 (md.results.TransientSolution(2).Thickness),... 80 (md.results.TransientSolution(2).MaskIceLevelset),... 81 (md.results.TransientSolution(2).CalvingMeltingrate),... 82 (md.results.TransientSolution(10).Vx),... 83 (md.results.TransientSolution(10).Vy),... 84 (md.results.TransientSolution(10).Vel),... 85 (md.results.TransientSolution(10).Pressure),... 86 (md.results.TransientSolution(10).Base),... 87 (md.results.TransientSolution(10).Surface),... 88 (md.results.TransientSolution(10).Thickness),... 89 (md.results.TransientSolution(10).MaskIceLevelset),... 90 (md.results.TransientSolution(10).CalvingMeltingrate),... 91 }; 102 (md.results.TransientSolution(1).BasalforcingsFloatingiceMeltingRate),... 103 (md.results.TransientSolution(5).Vx),... 104 (md.results.TransientSolution(5).Vy),... 105 (md.results.TransientSolution(5).Vel),... 106 (md.results.TransientSolution(5).Thickness),... 107 (md.results.TransientSolution(5).MaskIceLevelset),... 108 (md.results.TransientSolution(5).CalvingCalvingrate),... 109 (md.results.TransientSolution(5).CalvingMeltingrate),... 110 (md.results.TransientSolution(5).BasalforcingsFloatingiceMeltingRate),... 111 (md.results.TransientSolution(10).Vx),... 112 (md.results.TransientSolution(10).Vy),... 113 (md.results.TransientSolution(10).Vel),... 114 (md.results.TransientSolution(10).Thickness),... 115 (md.results.TransientSolution(10).MaskIceLevelset),... 116 (md.results.TransientSolution(10).CalvingCalvingrate),... 117 (md.results.TransientSolution(10).CalvingMeltingrate),... 118 (md.results.TransientSolution(10).BasalforcingsFloatingiceMeltingRate),... 119 }; -
../trunk-jpl/test/NightlyRun/test543.py
1 #Test Name: PigTran Frontalforcingsrignotautoregression1 #Test Name: PigTranStochasticforcing 2 2 import numpy as np 3 3 4 from calvingvonmises import *5 4 from frontalforcingsrignotautoregression import * 6 5 from socket import gethostname 7 6 from model import * … … 18 17 md = setflowequation(md, 'SSA', 'all') 19 18 md.timestepping.start_time = 0 20 19 md.timestepping.time_step = 1 21 md.timestepping.final_time = 2520 md.timestepping.final_time = 10 22 21 23 # Separate domain in 2 basins24 idb asin= np.zeros((md.mesh.numberofelements,))22 # Basin separation TF 23 idb_tf = np.zeros((md.mesh.numberofelements,)) 25 24 iid1 = np.where(md.mesh.x <= -1.6e6)[0] 26 25 for ii in range(md.mesh.numberofelements): 27 26 for vertex in range(3): 28 27 if md.mesh.elements[ii][vertex] - 1 in iid1: # one vertex in basin 1; NOTE: offset because of 1-based vertex indexing 29 idb asin[ii] = 128 idb_tf[ii] = 1 30 29 if idbasin[ii] == 0: # no vertex was found in basin 1 31 30 for vertex in range(3): 32 idbasin[ii] = 2 31 idb_tf[ii] = 2 32 # Basin separation default 33 idb_df = np.zeros((md.mesh.numberofelements,)) 34 iid1 = np.where(md.mesh.x <= -1.62e6)[0] 35 for ii in range(md.mesh.numberofelements): 36 for vertex in range(3): 37 if md.mesh.elements[ii][vertex] - 1 in iid1: # one vertex in basin 1; NOTE: offset because of 1-based vertex indexing 38 idb_df[ii] = 1 39 if idbasin[ii] == 0: # no vertex was found in basin 1 40 for vertex in range(3): 41 idb_df[ii] = 2 42 #Dimensionalities 43 nb_tf = 2 44 nb_clv = 2 45 nb_flmlt = 2 33 46 34 # Calving and frontalforcingsparameters47 # Calving parameters 35 48 md.mask.ice_levelset = 1e4 * (md.mask.ice_levelset + 0.5) 36 md.calving = calvingvonmises() 49 md.calving.calvingrate = 0.3*np.ones((md.mesh.numberofvertices,)) 50 md.levelset.spclevelset = np.full((md.mesh.numberofvertices,), np.nan) 51 md.levelset.migration_max = 10.0 52 #Frontal forcing parameters 37 53 md.frontalforcings = frontalforcingsrignotautoregression() 38 md.frontalforcings.num_basins = 239 md.frontalforcings.basin_id = idb asin- 1 # NOTE: offset because of 1-based vertex indexing54 md.frontalforcings.num_basins = nb_tf 55 md.frontalforcings.basin_id = idb_tf - 1 # NOTE: offset because of 1-based vertex indexing 40 56 md.frontalforcings.subglacial_discharge = 0.1 * np.ones((md.mesh.numberofvertices,)) 41 57 md.frontalforcings.beta0 = np.array([[0.05, 0.01]]) # intercept values of TF in basins [C] 42 58 md.frontalforcings.beta1 = np.array([[0.001, 0.0001]]) # trend values of TF in basins [C/yr] … … 44 60 md.frontalforcings.ar_order = 4 45 61 md.frontalforcings.ar_timestep = 2 # timestep of the autoregressive model [yr] 46 62 md.frontalforcings.phi = np.array([[0.1, -0.1, 0.01, -0.01], [0.2, -0.2, 0.1, 0.0]]) # autoregressive parameters 63 #Floating Ice Melt parameters 64 md.basalforcings.floatingice_melting_rate = 0.1*np.ones((md.mesh.numberofvertices,)) 47 65 66 #Covariance matrix 67 covtf = 1e-4*np.identity(nb_tf) 68 covclv = 1e-1*np.identity(nb_clv) 69 covclv[0,0] = 1/10*covclv[0,0] 70 covflmlt = 0.05*np.identity(nb_flmlt) 71 covglob = np.zeros([6,6]) 72 covglob[0:2,0:2] = covtf 73 covglob[2:4,2:4] = covclv 74 covglob[4:6,4:6] = covflmlt 75 48 76 # Stochastic forcing 49 77 md.stochasticforcing.isstochasticforcing = 1 50 md.stochasticforcing.fields = ['FrontalForcingsRignotAutoregression'] 78 md.stochasticforcing.fields = ['FrontalForcingsRignotAutoregression','DefaultCalving','FloatingMeltRate'] 79 md.stochasticforcing.defauldimension = 2 80 md.stochasticforcing.default_id = idb_df 81 md.stochasticforcing.covariance = covglob # global covariance among- and between-fields 82 md.stochasticforcing.randomflag = 0 # determines true/false randomness 51 83 52 md. stochasticforcing.covariance = 1e-4 * np.array([[1.5, 0.5], [0.5, 0.4]]) # global covariance among- and between-fields53 md. stochasticforcing.randomflag = 0 # determines true/false randomness84 md.transient.ismovingfront = 1 85 md.transient.isgroundingline = 1 54 86 55 56 md.transient.ismovingfront = 1 57 md.levelset.spclevelset = np.full((md.mesh.numberofvertices,), np.nan) 58 md.levelset.migration_max = 1e10 59 60 md.transient.requested_outputs = ['default', 'CalvingMeltingrate'] 87 md.transient.requested_outputs = ['default', 'CalvingCalvingrate','CalvingMeltingrate','BasalforcingsFloatingiceMeltingRate'] 61 88 md.cluster = generic('name', gethostname(), 'np', 2) 62 89 md = solve(md, 'Transient') 63 90 64 91 # Fields and tolerances to track changes 65 92 field_names = [ 66 'Vx1' , 'Vy1', 'Vel1', 'Pressure1', 'Bed1', 'Surface1', 'Thickness1', 'MaskIceLevelset1', 'CalvingMetlingRate1',67 'Vx2' , 'Vy2', 'Vel2', 'Pressure2', 'Bed2', 'Surface2', 'Thickness2', 'MaskIceLevelset2', 'CalvingMetlingRate2',68 'Vx10', 'Vy10', 'Vel10', 'Pressure10', 'Bed10', 'Surface10', 'Thickness10', 'MaskIceLevelset10', 'CalvingMetlingRate10'93 'Vx1' ,'Vy1' ,'Vel1' ,'Thickness1' ,'MaskIceLevelset1' ,'CalvingCalvingrate1' ,'CalvingMeltingrate1' ,'BasalforcingsFloatingiceMeltingRate1', 94 'Vx2' ,'Vy2' ,'Vel2' ,'Thickness2' ,'MaskIceLevelset2' ,'CalvingCalvingrate1' ,'CalvingMeltingrate2' ,'BasalforcingsFloatingiceMeltingRate2', 95 'Vx10','Vy10','Vel10','Thickness10','MaskIceLevelset10','CalvingCalvingrate10','CalvingMeltingrate10','BasalforcingsFloatingiceMeltingRate10' 69 96 ] 97 70 98 field_tolerances = [ 71 1e-11, 2e-11, 2e-11, 1e-12, 2e-11, 6e-12, 9e-12, 1e-12, 1e-9,72 2e-11, 1e-11, 1e-11, 9e-12, 2e-1, 2e-11, 2e-11, 1e-11, 1e-9,73 2e-10, 1e-10, 1e-10, 1e-10, 1e-10, 1e-10, 1e-10, 1e-9, 1e-999 1e-11,2e-11,2e-11,1e-11,1e-9,1e-10,1e-10,1e-10, 100 2e-11,1e-11,1e-11,9e-11,2e-9,1e-10,1e-10,1e-10, 101 2e-10,1e-10,1e-10,1e-10,5e-9,1e-10,1e-10,1e-10 74 102 ] 75 103 field_values = [ 76 104 md.results.TransientSolution[0].Vx, 77 105 md.results.TransientSolution[0].Vy, 78 106 md.results.TransientSolution[0].Vel, 79 md.results.TransientSolution[0].Pressure,80 md.results.TransientSolution[0].Base,81 md.results.TransientSolution[0].Surface,82 107 md.results.TransientSolution[0].Thickness, 83 108 md.results.TransientSolution[0].MaskIceLevelset, 109 md.results.TransientSolution[0].CalvingCalvingrate, 84 110 md.results.TransientSolution[0].CalvingMeltingrate, 85 md.results.TransientSolution[ 1].Vx,86 md.results.TransientSolution[ 1].Vy,87 md.results.TransientSolution[ 1].Vel,88 md.results.TransientSolution[ 1].Pressure,89 md.results.TransientSolution[ 1].Base,90 md.results.TransientSolution[ 1].Surface,91 md.results.TransientSolution[ 1].Thickness,92 md.results.TransientSolution[ 1].MaskIceLevelset,93 md.results.TransientSolution[ 1].CalvingMeltingrate,111 md.results.TransientSolution[0].BasalforcingsFloatingiceMeltingRate, 112 md.results.TransientSolution[4].Vx, 113 md.results.TransientSolution[4].Vy, 114 md.results.TransientSolution[4].Vel, 115 md.results.TransientSolution[4].Thickness, 116 md.results.TransientSolution[4].MaskIceLevelset, 117 md.results.TransientSolution[4].CalvingCalvingrate, 118 md.results.TransientSolution[4].CalvingMeltingrate, 119 md.results.TransientSolution[4].BasalforcingsFloatingiceMeltingRate, 94 120 md.results.TransientSolution[9].Vx, 95 121 md.results.TransientSolution[9].Vy, 96 122 md.results.TransientSolution[9].Vel, 97 md.results.TransientSolution[9].Pressure,98 md.results.TransientSolution[9].Base,99 md.results.TransientSolution[9].Surface,100 123 md.results.TransientSolution[9].Thickness, 101 124 md.results.TransientSolution[9].MaskIceLevelset, 102 md.results.TransientSolution[9].CalvingMeltingrate 125 md.results.TransientSolution[9].CalvingCalvingrate, 126 md.results.TransientSolution[9].CalvingMeltingrate, 127 md.results.TransientSolution[9].BasalforcingsFloatingiceMeltingRate 103 128 ]
Note:
See TracBrowser
for help on using the repository browser.