[4949] | 1 | %Start defining model parameters here
|
---|
| 2 |
|
---|
[4953] | 3 | %Geometry and observation
|
---|
[21341] | 4 | x = archread('../Data/Pig.arch','x');
|
---|
| 5 | y = archread('../Data/Pig.arch','y');
|
---|
| 6 | vx_obs = archread('../Data/Pig.arch','vx_obs');
|
---|
| 7 | vy_obs = archread('../Data/Pig.arch','vy_obs');
|
---|
| 8 | index = archread('../Data/Pig.arch','index');
|
---|
| 9 | surface = archread('../Data/Pig.arch','surface');
|
---|
| 10 | thickness = archread('../Data/Pig.arch','thickness');
|
---|
[24313] | 11 | bed = archread('../Data/Pig.arch','bed');
|
---|
[21341] | 12 |
|
---|
| 13 | x = x{1};
|
---|
| 14 | y = y{1};
|
---|
| 15 | vx_obs = vx_obs{1};
|
---|
| 16 | vy_obs = vy_obs{1};
|
---|
| 17 | index = index{1};
|
---|
| 18 | surface = surface{1};
|
---|
| 19 | thickness = thickness{1};
|
---|
[24313] | 20 | bed = bed{1};
|
---|
[21341] | 21 |
|
---|
[9734] | 22 | md.inversion.vx_obs =InterpFromMeshToMesh2d(index,x,y,vx_obs,md.mesh.x,md.mesh.y);
|
---|
| 23 | md.inversion.vy_obs =InterpFromMeshToMesh2d(index,x,y,vy_obs,md.mesh.x,md.mesh.y);
|
---|
| 24 | md.geometry.surface =InterpFromMeshToMesh2d(index,x,y,surface,md.mesh.x,md.mesh.y);
|
---|
| 25 | md.geometry.thickness=InterpFromMeshToMesh2d(index,x,y,thickness,md.mesh.x,md.mesh.y);
|
---|
[17806] | 26 | md.geometry.base=md.geometry.surface-md.geometry.thickness;
|
---|
[24313] | 27 | md.geometry.bed =md.geometry.base;
|
---|
[25836] | 28 | pos = find(md.mask.ocean_levelset<0);
|
---|
[24313] | 29 | md.geometry.bed(pos) =InterpFromMeshToMesh2d(index,x,y,bed,md.mesh.x(pos),md.mesh.y(pos));
|
---|
[4953] | 30 | clear surface thickness vx_obs vy_obs x y index;
|
---|
[9684] | 31 | md.initialization.vx=md.inversion.vx_obs;
|
---|
| 32 | md.initialization.vy=md.inversion.vy_obs;
|
---|
[9725] | 33 | md.initialization.vz=zeros(md.mesh.numberofvertices,1);
|
---|
| 34 | md.initialization.pressure=zeros(md.mesh.numberofvertices,1);
|
---|
[4949] | 35 |
|
---|
[4953] | 36 | %Materials
|
---|
[9725] | 37 | md.initialization.temperature=(273-20)*ones(md.mesh.numberofvertices,1);
|
---|
[9684] | 38 | md.materials.rheology_B=paterson(md.initialization.temperature);
|
---|
[9725] | 39 | md.materials.rheology_n=3*ones(md.mesh.numberofelements,1);
|
---|
[9684] | 40 | md.initialization.temperature=md.initialization.temperature;
|
---|
[4953] | 41 |
|
---|
| 42 | %Friction
|
---|
[9725] | 43 | md.friction.coefficient=50*ones(md.mesh.numberofvertices,1);
|
---|
[25836] | 44 | md.friction.coefficient(find(md.mask.ocean_levelset<0.))=0.;
|
---|
[9725] | 45 | md.friction.p=ones(md.mesh.numberofelements,1);
|
---|
| 46 | md.friction.q=ones(md.mesh.numberofelements,1);
|
---|
[4953] | 47 |
|
---|
| 48 | %Numerical parameters
|
---|
[16137] | 49 | md.masstransport.stabilization=1;
|
---|
[6307] | 50 | md.verbose=verbose(0);
|
---|
[9702] | 51 | md.settings.waitonlock=30;
|
---|
[9628] | 52 | md.timestepping.time_step=1;
|
---|
| 53 | md.timestepping.final_time=2;
|
---|
[16137] | 54 | md.stressbalance.restol=0.05;
|
---|
| 55 | md.stressbalance.reltol=1;
|
---|
[9677] | 56 | md.steadystate.reltol=1;
|
---|
[16137] | 57 | md.stressbalance.abstol=NaN;
|
---|
[23189] | 58 | md.groundingline.migration='None';
|
---|
[4949] | 59 |
|
---|
[4953] | 60 | %Boundary conditions:
|
---|
[4949] | 61 | md=SetMarineIceSheetBC(md);
|
---|
[4999] | 62 |
|
---|
| 63 | %Change name so that no test have the same name
|
---|
| 64 | A=dbstack;
|
---|
[9625] | 65 | if (length(A)>2), md.miscellaneous.name=A(3).file(1:end-2); end
|
---|