Changeset 20741

Timestamp:

06/15/16 17:41:06 (9 years ago)

Author:

schlegel

Message:

CHG: update for workshop

Location:

issm/trunk-jpl/examples/Greenland

Files:

• CheatSheet.m (modified) (1 diff)
• Greenland.par (modified) (1 diff)
• runme.m (modified) (6 diffs)

issm/trunk-jpl/examples/Greenland/CheatSheet.m

@@ +1 @@
+if any(steps==7)
+        disp('   Step 7: Historical Relaxation run');
+        md = loadmodel('./Models/Greenland.Control_drag');
+
+        load smbbox
+
+        %convert mesh x,y into the Box projection
+        [md.mesh.lat,md.mesh.long]  = xy2ll(md.mesh.x,md.mesh.y,+1,39,71);
+        [xi,yi]= ll2xy(md.mesh.lat,md.mesh.long,+1,45,70);
+
+        %Interpolate and set surface mass balance
+        index = BamgTriangulate(x1(:),y1(:));
+        smb_mo = InterpFromMeshToMesh2d(index,x1(:),y1(:),smbmean(:),xi,yi);
+        smb = smb_mo*12/1000*md.materials.rho_freshwater/md.materials.rho_ice;
+        md.smb.mass_balance = [smb;1 ];
+
+        %Set transient options, run for 20 years, saving every 5 timesteps
+        md.timestepping.time_step=0.2;
+        md.timestepping.final_time=200;
+        md.settings.output_frequency=5;
+
+        %Additional options
+        md.inversion.iscontrol=0;
+        md.transient.requested_outputs={'IceVolume','TotalSmb', ...
+                'SmbMassBalance'};
+        md.verbose=verbose('solution',true,'module',true);
+
+        %Go solve
+        md.cluster=generic('name',oshostname,'np',2);
+        md=solve(md,TransientSolutionEnum);
+
+        save ./Models/Greenland.HistoricTransient_200yr md;
+end
+
 if any(steps==8)
         disp('   Step 8: Plotting exercise');
 
         %Load historic transient model
-        md = loadmodel('./Models/Greenland.HistoricTransient');
+        md = loadmodel('./Models/Greenland.HistoricTransient_200yr');
 
         %Create Line Plots of relaxation run.  Create a figure.
         figure
 
-        %Save surface mass balance, by looping through 200 years, or 1000 steps
-        surfmb=[]; for i=1:1000; surfmb=[surfmb ...
+        %Save surface mass balance, by looping through 200 years (1000 steps)
+        % Note, the first output will always contain output from time step 1
+        surfmb=[]; for i=2:201; surfmb=[surfmb ...
                 md.results.TransientSolution(i).SmbMassBalance]; end
 
         %Plot surface mass balance time series in first subplot
-        subplot(3,1,1); plot([0.2:0.2:200],mean(surfmb));
+        subplot(3,1,1); plot([1:200],mean(surfmb));
 
         %Title this plot Mean surface mass balance
         title('Mean Surface mass balance');
 
-        %Save velocity by looping through 1000 steps
-        vel=[]; for i=1:1000; vel=[vel md.results.TransientSolution(i).Vel]; end
+        %Save velocity by looping through 200 years
+        vel=[]; for i=2:201; vel=[vel md.results.TransientSolution(i).Vel]; end
 
         %Plot velocity time series in second subplot
-        subplot(3,1,2); plot([0.2:0.2:200],mean(vel));
+        subplot(3,1,2); plot([1:200],mean(vel));
 
         %Title this plot Mean Velocity
         title('Mean Velocity');
 
-        %Save Ice Volume by looping through 1000 steps
-        volume=[]; for i=1:1000; volume=[volume md.results.TransientSolution(i).IceVolume]; end
+        %Save Ice Volume by looping through 200 years
+        volume=[]; for i=2:201; volume=[volume md.results.TransientSolution(i).IceVolume]; end
 
         %Plot volume time series in third subplot
-        subplot(3,1,3); plot([0.2:0.2:200],volume);
+        subplot(3,1,3); plot([1:200],volume);
 
         %Title this plot Mean Velocity and add an x label of years

issm/trunk-jpl/examples/Greenland/Greenland.par

@@ -69 +69 @@
 md.initialization.pressure=md.materials.rho_ice*md.constants.g*md.geometry.thickness;
 
-disp('   Single point constraint for continental model');
-%Set at least one vertex to velocity 0 so as to not get a singular problem (point on the wet peninsula)
+disp('   Single point constraints');
+%Initialize single point constraint arrarys
 md.stressbalance.referential=NaN*ones(md.mesh.numberofvertices,6);
 md.stressbalance.spcvx = NaN*ones(md.mesh.numberofvertices,1);
 md.stressbalance.spcvy = NaN*ones(md.mesh.numberofvertices,1);
 md.stressbalance.spcvz = NaN*ones(md.mesh.numberofvertices,1);
-location = 1.0e+06 *[.32011 -2.2039];
-[dist pos]=min(sqrt((md.mesh.x - location(1)).^2 + (md.mesh.y - location(2)).^2));
-md.stressbalance.spcvx(pos) = 0;
-md.stressbalance.spcvy(pos) = 0;
-md.stressbalance.spcvz(pos) = 0;
 

issm/trunk-jpl/examples/Greenland/runme.m

@@ -117 +117 @@
 
         %Planview plots
-        plotmodel(md,'data',md.results.TransientSolution(end).Vel,'caxis',[1e-1 6000],'log', 10, 'title', 'Velocity (m/y)','gridded',11);
+        plotmodel(md,'data',md.results.TransientSolution(end).Vel,'log#1', 10, ...
+                'caxis#1', [1e-1 1e4], ...
+                'title#1', 'Velocity (m/y)',...
+                'data',md.results.TransientSolution(1).SmbMassBalance,...
+                'title#2', 'Surface Mass Balance (m/y)',...
+                'data',md.results.TransientSolution(end).Thickness,...
+                'title', 'Thickness (m)',...
+                'data',md.results.TransientSolution(end).Surface,...
+                'title', 'Surface (m)');
 
         %Line Plots
@@ -161 +169 @@
 
         %plot a time series of total SMB
-        plot(t,totalsmb/1e9); title('Total Surface mass balance, Gt'); xlabel('year'); ylabel('Gt/yr');
+        figure; plot(t,totalsmb/1e9); title('Total Surface mass balance, Gt'); xlabel('year'); ylabel('Gt/yr');
 
         clear smbbox
@@ -210 +218 @@
 
 
-        %Save surface mass balance, by looping through 200 years, or 1000 steps
+        %Save surface mass balance, by looping through 200 years (1000 steps)
+        % Note, the first output will always contain output from time step 1
 
 
@@ -219 +228 @@
 
 
-        %Save velocity by looping through 1000 steps
+        %Save velocity by looping through 200 years
 
 
@@ -228 +237 @@
 
 
-        %Save Ice Volume by looping through 1000 steps
+        %Save Ice Volume by looping through 200 years 
 
 
@@ -240 +249 @@
 if any(steps==9)
         disp('   Step 9: Box Transient run');
-        md = loadmodel('./Models/Greenland.HistoricTransient');
+        md = loadmodel('./Models/Greenland.HistoricTransient_200yr');
 
         %load past transient results