1 | %Start defining model parameters here
|
---|
2 |
|
---|
3 | %Geometry and observation
|
---|
4 | x = transpose(ncread('../Data/Pig.nc','x'));
|
---|
5 | y = transpose(ncread('../Data/Pig.nc','y'));
|
---|
6 | vx_obs = transpose(ncread('../Data/Pig.nc','vx_obs'));
|
---|
7 | vy_obs = transpose(ncread('../Data/Pig.nc','vy_obs'));
|
---|
8 | index = transpose(ncread('../Data/Pig.nc','index'));
|
---|
9 | surface = transpose(ncread('../Data/Pig.nc','surface'));
|
---|
10 | thickness = transpose(ncread('../Data/Pig.nc','thickness'));
|
---|
11 | md.inversion.vx_obs =InterpFromMeshToMesh2d(index,x,y,vx_obs,md.mesh.x,md.mesh.y);
|
---|
12 | md.inversion.vy_obs =InterpFromMeshToMesh2d(index,x,y,vy_obs,md.mesh.x,md.mesh.y);
|
---|
13 | md.geometry.surface =InterpFromMeshToMesh2d(index,x,y,surface,md.mesh.x,md.mesh.y);
|
---|
14 | md.geometry.thickness=InterpFromMeshToMesh2d(index,x,y,thickness,md.mesh.x,md.mesh.y);
|
---|
15 | md.geometry.base=md.geometry.surface-md.geometry.thickness;
|
---|
16 | clear surface thickness vx_obs vy_obs x y index;
|
---|
17 | md.initialization.vx=md.inversion.vx_obs;
|
---|
18 | md.initialization.vy=md.inversion.vy_obs;
|
---|
19 | md.initialization.vz=zeros(md.mesh.numberofvertices,1);
|
---|
20 | md.initialization.pressure=zeros(md.mesh.numberofvertices,1);
|
---|
21 |
|
---|
22 | %Materials
|
---|
23 | md.initialization.temperature=(273-20)*ones(md.mesh.numberofvertices,1);
|
---|
24 | md.materials.rheology_B=paterson(md.initialization.temperature);
|
---|
25 | md.materials.rheology_n=3*ones(md.mesh.numberofelements,1);
|
---|
26 | md.initialization.temperature=md.initialization.temperature;
|
---|
27 |
|
---|
28 | %Friction
|
---|
29 | md.friction.coefficient=50*ones(md.mesh.numberofvertices,1);
|
---|
30 | md.friction.coefficient(find(md.mask.groundedice_levelset<0.))=0.;
|
---|
31 | md.friction.p=ones(md.mesh.numberofelements,1);
|
---|
32 | md.friction.q=ones(md.mesh.numberofelements,1);
|
---|
33 |
|
---|
34 | %Numerical parameters
|
---|
35 | md.stressbalance.viscosity_overshoot=0.3;
|
---|
36 | md.masstransport.stabilization=1;
|
---|
37 | md.verbose=verbose(0);
|
---|
38 | md.settings.waitonlock=30;
|
---|
39 | md.timestepping.time_step=1;
|
---|
40 | md.timestepping.final_time=2;
|
---|
41 | md.stressbalance.restol=0.05;
|
---|
42 | md.stressbalance.reltol=1;
|
---|
43 | md.steadystate.reltol=1;
|
---|
44 | md.stressbalance.abstol=NaN;
|
---|
45 |
|
---|
46 | %Boundary conditions:
|
---|
47 | md=SetMarineIceSheetBC(md);
|
---|
48 |
|
---|
49 | %Change name so that no test have the same name
|
---|
50 | A=dbstack;
|
---|
51 | if (length(A)>2), md.miscellaneous.name=A(3).file(1:end-2); end
|
---|