Changeset 20746
- Timestamp:
- 06/15/16 19:08:16 (9 years ago)
- Location:
- issm/trunk-jpl/examples/IceBridge
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
issm/trunk-jpl/examples/IceBridge/Greenland.par
r20743 r20746 26 26 md.geometry.thickness(pos0)=1; 27 27 md.geometry.surface=md.geometry.thickness+md.geometry.base; 28 29 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 30 %Reading IceBridge data for Jakobshavn 31 disp(' reading IceBridge Jakobshavn bedrock'); 32 fid = fopen('../Data/Jakobshavn_2008_2011_Composite_XYZGrid.txt'); 33 titles = fgets(fid); 34 data = fscanf(fid,'%g,%g,%g,%g,%g',[5 266400])'; 35 fclose(fid); 36 37 [xi,yi]= ll2xy(md.mesh.lat,md.mesh.long,+1,45,70); 38 bed = flipud(reshape(data(:,5),[360 740])); bed(find(bed==-9999))=NaN; 39 bedy = flipud(reshape(data(:,1),[360 740])); 40 bedx = flipud(reshape(data(:,2),[360 740])); 41 42 %Insert Icebridge bed and recalculate thickness 43 bed_jks=InterpFromGridToMesh(bedx(1,:)',bedy(:,1),bed,xi,yi,NaN); 44 in=ContourToMesh(md.mesh.elements,md.mesh.x,md.mesh.y,... 45 'Jak_grounded.exp','node',1); 46 bed_jks(~in)=NaN; 47 pos=find(~isnan(bed_jks)); 48 md.geometry.base(pos)=bed_jks(pos); 49 md.geometry.thickness=md.geometry.surface-md.geometry.base; 50 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 28 51 29 52 disp(' Interpolating velocities '); … … 58 81 disp(' Set other boundary conditions'); 59 82 md.mask.ice_levelset(md.mesh.vertexonboundary==1)=0; 83 md.basalforcings.groundedice_melting_rate = zeros(md.mesh.numberofvertices,1); 60 84 md.basalforcings.floatingice_melting_rate = zeros(md.mesh.numberofvertices,1); 61 md.basalforcings.groundedice_melting_rate = zeros(md.mesh.numberofvertices,1);62 85 md.thermal.spctemperature = [md.initialization.temperature;1]; %impose observed temperature on surface 63 86 md.masstransport.spcthickness = NaN*ones(md.mesh.numberofvertices,1); … … 70 93 71 94 disp(' Single point constraints'); 72 %Initialize single point constraint arra rys95 %Initialize single point constraint arrays 73 96 md.stressbalance.referential=NaN*ones(md.mesh.numberofvertices,6); 74 97 md.stressbalance.spcvx = NaN*ones(md.mesh.numberofvertices,1); 75 98 md.stressbalance.spcvy = NaN*ones(md.mesh.numberofvertices,1); 76 99 md.stressbalance.spcvz = NaN*ones(md.mesh.numberofvertices,1); 77 -
issm/trunk-jpl/examples/IceBridge/Greenland_solution.par
r20743 r20746 1 %Name and hemisphere1 %Name and Coordinate system 2 2 md.miscellaneous.name='SeaRISEgreenland'; 3 3 md.mesh.epsg=3413; 4 4 5 5 disp(' Loading SeaRISE data from NetCDF'); … … 31 31 disp(' reading IceBridge Jakobshavn bedrock'); 32 32 fid = fopen('../Data/Jakobshavn_2008_2011_Composite_XYZGrid.txt'); 33 titles = fgets(fid); data = fscanf(fid,'%g,%g,%g,%g,%g',[5 266400])'; 33 titles = fgets(fid); 34 data = fscanf(fid,'%g,%g,%g,%g,%g',[5 266400])'; 34 35 fclose(fid); 35 36 … … 44 45 surf_jks=InterpFromGridToMesh(bedx(1,:)',bedy(:,1),surf,xi,yi,NaN); 45 46 in=ContourToMesh(md.mesh.elements,md.mesh.x,md.mesh.y,... 46 'Jak_grounded.exp','node',1);47 'Jak_grounded.exp','node',1); 47 48 bed_jks(~in)=NaN; 48 49 surf_jks(~in)=NaN; … … 66 67 67 68 disp(' Interpolating surface mass balance'); 68 md.s urfaceforcings.mass_balance=InterpFromGridToMesh(x1,y1,smb,md.mesh.x,md.mesh.y,0);69 md.s urfaceforcings.mass_balance=md.surfaceforcings.mass_balance*md.materials.rho_water/md.materials.rho_ice;69 md.smb.mass_balance=InterpFromGridToMesh(x1,y1,smb,md.mesh.x,md.mesh.y,0); 70 md.smb.mass_balance=md.smb.mass_balance*md.materials.rho_water/md.materials.rho_ice; 70 71 71 72 disp(' Construct basal friction parameters'); … … 84 85 disp(' Set other boundary conditions'); 85 86 md.mask.ice_levelset(md.mesh.vertexonboundary==1)=0; 86 md.basalforcings.melting_rate = zeros(md.mesh.numberofvertices,1); 87 md.basalforcings.groundedice_melting_rate = zeros(md.mesh.numberofvertices,1); 88 md.basalforcings.floatingice_melting_rate = zeros(md.mesh.numberofvertices,1); 87 89 md.thermal.spctemperature = [md.initialization.temperature;1]; %impose observed temperature on surface 88 90 md.masstransport.spcthickness = NaN*ones(md.mesh.numberofvertices,1); … … 94 96 md.initialization.pressure=md.materials.rho_ice*md.constants.g*md.geometry.thickness; 95 97 96 disp(' Single point constraint for continental model');97 % Set at least one vertex to velocity 0 so as to not get a singular problem (point on the wet peninsula)98 disp(' Single point constraints'); 99 %Initialize single point constraint arrays 98 100 md.stressbalance.referential=NaN*ones(md.mesh.numberofvertices,6); 99 101 md.stressbalance.spcvx = NaN*ones(md.mesh.numberofvertices,1); 100 102 md.stressbalance.spcvy = NaN*ones(md.mesh.numberofvertices,1); 101 103 md.stressbalance.spcvz = NaN*ones(md.mesh.numberofvertices,1); 102 location = 1.0e+06 *[.32011 -2.2039];103 [dist pos]=min(sqrt((md.mesh.x - location(1)).^2 + (md.mesh.y - location(2)).^2));104 md.stressbalance.spcvx(pos) = 0;105 md.stressbalance.spcvy(pos) = 0;106 md.stressbalance.spcvz(pos) = 0;107
Note:
See TracChangeset
for help on using the changeset viewer.