1 | %Start defining model parameters here
|
---|
2 |
|
---|
3 | %Geometry
|
---|
4 | hmin=300;
|
---|
5 | hmax=1000;
|
---|
6 | ymin=min(md.y);
|
---|
7 | ymax=max(md.y);
|
---|
8 | md.geometry.thickness=hmax+(hmin-hmax)*(md.y-ymin)/(ymax-ymin);
|
---|
9 | md.geometry.bed=-md.materials.rho_ice/md.materials.rho_water*md.geometry.thickness;
|
---|
10 | bed_sheet=-md.materials.rho_ice/md.materials.rho_water*(hmax+(hmin-hmax)*(ymax/2-ymin)/(ymax-ymin));
|
---|
11 | pos=find(md.y<=ymax/2);
|
---|
12 | md.geometry.bed(pos)=bed_sheet;
|
---|
13 | md.geometry.surface=md.geometry.bed+md.geometry.thickness;
|
---|
14 |
|
---|
15 | %Initial velocity
|
---|
16 | load('../Data/SquareSheetShelf.data','-mat');
|
---|
17 | md.initialization.vx=InterpFromMeshToMesh2d(index,x,y,vx,md.x,md.y);
|
---|
18 | md.initialization.vy=InterpFromMeshToMesh2d(index,x,y,vy,md.x,md.y);
|
---|
19 | clear vx vy x y index;
|
---|
20 | md.initialization.vz=zeros(md.mesh.numberofvertices,1);
|
---|
21 | md.initialization.pressure=zeros(md.mesh.numberofvertices,1);
|
---|
22 |
|
---|
23 | %Materials
|
---|
24 | md.initialization.temperature=(273-20)*ones(md.mesh.numberofvertices,1);
|
---|
25 | md.materials.rheology_B=paterson(md.initialization.temperature);
|
---|
26 | md.materials.rheology_n=3*ones(md.mesh.numberofelements,1);
|
---|
27 |
|
---|
28 | %Accumulation and melting
|
---|
29 | md.surfaceforcings.mass_balance=10*ones(md.mesh.numberofvertices,1);
|
---|
30 | md.basalforcings.melting_rate=5*ones(md.mesh.numberofvertices,1);
|
---|
31 |
|
---|
32 | %Friction
|
---|
33 | pos=find(md.mask.elementonfloatingice);
|
---|
34 | md.friction.coefficient=20*ones(md.mesh.numberofvertices,1);
|
---|
35 | md.friction.coefficient(md.mesh.elements(pos,:))=0;
|
---|
36 | md.friction.p=ones(md.mesh.numberofelements,1);
|
---|
37 | md.friction.q=ones(md.mesh.numberofelements,1);
|
---|
38 |
|
---|
39 | %Numerical parameters
|
---|
40 | md.diagnostic.viscosity_overshoot=0.0;
|
---|
41 | md.prognostic.stabilization=1;
|
---|
42 | md.thermal.stabilization=1;
|
---|
43 | md.verbose=verbose(0);
|
---|
44 | md.settings.waitonlock=30;
|
---|
45 | md.diagnostic.restol=0.05;
|
---|
46 | md.steadystate.reltol=0.05;
|
---|
47 | md.diagnostic.reltol=0.05;
|
---|
48 | md.diagnostic.abstol=NaN;
|
---|
49 | md.timestepping.time_step=1;
|
---|
50 | md.timestepping.final_time=3;
|
---|
51 |
|
---|
52 | %Deal with boundary conditions:
|
---|
53 | md=SetMarineIceSheetBC(md,'../Exp/SquareFront.exp');
|
---|
54 |
|
---|
55 | %Change name so that no test have the same name
|
---|
56 | A=dbstack;
|
---|
57 | if (length(A)>2), md.miscellaneous.name=A(3).file(1:end-2); end
|
---|