<p><b>mhoffman@lanl.gov</b> 2013-01-30 20:40:57 -0700 (Wed, 30 Jan 2013)</p><p>BRANCH COMMIT - land ice<br>
<br>
Some first modifications necessary for modeling ice shelves.<br>
<br>
1. Added a marine basal mass balance field. This is input as a time series (but can only have 1 time level) and is defined across the entire domain. Each time step the mask is used to 0 this field for grounded ice cells. Therefore this field can be defined for areas that will unground in the future. For now this is simply an input field, but eventually it can come from a coupler or ocean model. Note: grounded ice basal mass balance is not yet defined because a column diffusion temperature solve does not yet exist.<br>
<br>
2. Set beta to 0 for floating ice. LifeV may already do this, but doing it here for completeness.<br>
<br>
3. Adjusting Registry to output lower and upper surface which are useful for visualizing floating ice.<br>
<br>
(also added messages to stderr when trying to run a lifev velocity solver without a lifev build. such a message would have saved me some time recently.)<br>
</p><hr noshade><pre><font color="gray">Modified: branches/land_ice_projects/implement_core/src/core_land_ice/Registry
===================================================================
--- branches/land_ice_projects/implement_core/src/core_land_ice/Registry        2013-01-30 22:34:57 UTC (rev 2399)
+++ branches/land_ice_projects/implement_core/src/core_land_ice/Registry        2013-01-31 03:40:57 UTC (rev 2400)
@@ -87,8 +87,8 @@
dim nBetaTimeSlices nBetaTimeSlices
dim nSfcAirTempTimeSlices nSfcAirTempTimeSlices
dim nBasalHeatFluxTimeSlices nBasalHeatFluxTimeSlices
+dim nMarineBasalMassBalTimeSlices nMarineBasalMassBalTimeSlices
-
% VARIABLES =====================
%
% var persistence type name_in_file ( dims ) time_levs iro- name_in_code struct super-array array_class
@@ -191,14 +191,15 @@
var persistent real basalHeatFlux ( nCells ) 0 - basalHeatFlux mesh - -
var persistent real sfcAirTemp ( nCells ) 0 - sfcAirTemp mesh - -
var persistent real sfcMassBal ( nCells ) 0 - sfcMassBal mesh - -
+var persistent real marineBasalMassBal ( nCells ) 0 - marineBasalMassBal mesh - -
% Time-varying forcing for B.C.
var persistent real betaTimeSeries ( nBetaTimeSlices nCells ) 0 ir betaTimeSeries mesh - -
var persistent real basalHeatFluxTimeSeries ( nBasalHeatFluxTimeSlices nCells ) 0 ir basalHeatFluxTimeSeries mesh - -
var persistent real sfcAirTempTimeSeries ( nSfcAirTempTimeSlices nCells ) 0 ir sfcAirTempTimeSeries mesh - -
var persistent real sfcMassBalTimeSeries ( nSfcMassBalTimeSlices nCells ) 0 ir sfcMassBalTimeSeries mesh - -
+var persistent real marineBasalMassBalTimeSeries ( nMarineBasalMassBalTimeSlices nCells ) 0 ir marineBasalMassBalTimeSeries mesh - -
-
% STATE VARIABLES =====================
% Prognostic variables: read from input, saved in restart, and written to output
@@ -222,8 +223,8 @@
var persistent real uReconstructMeridional ( nVertLevels nCells Time ) 2 - uReconstructMeridional state - -
var persistent real layerThickness ( nVertLevels nCells Time ) 2 - layerThickness state - -
var persistent real layerThicknessEdge ( nVertLevels nEdges Time ) 2 - layerThicknessEdge state - -
-var persistent real lowerSurface ( nCells Time ) 2 - lowerSurface state - -
-var persistent real upperSurface ( nCells Time ) 2 - upperSurface state - -
+var persistent real lowerSurface ( nCells Time ) 2 o lowerSurface state - -
+var persistent real upperSurface ( nCells Time ) 2 o upperSurface state - -
var persistent real temperatureUpperBC ( nCells Time ) 2 - temperatureUpperBC state - -
var persistent real temperatureLowerBC ( nCells Time ) 2 - temperatureLowerBC state - -
% Other diagnostic variables: neither read nor written to any files
Modified: branches/land_ice_projects/implement_core/src/core_land_ice/mpas_land_ice_annual_forcing.F
===================================================================
--- branches/land_ice_projects/implement_core/src/core_land_ice/mpas_land_ice_annual_forcing.F        2013-01-30 22:34:57 UTC (rev 2399)
+++ branches/land_ice_projects/implement_core/src/core_land_ice/mpas_land_ice_annual_forcing.F        2013-01-31 03:40:57 UTC (rev 2400)
@@ -94,8 +94,8 @@
! local variables
!
!-----------------------------------------------------------------
- real (kind=RKIND), dimension(:,:), pointer :: sfcMassBalTimeSeries, sfcAirTempTimeSeries, betaTimeSeries, basalHeatFluxTimeSeries
- real (kind=RKIND), dimension(:), pointer :: sfcMassBal, sfcAirTemp, beta, basalHeatFlux
+ real (kind=RKIND), dimension(:,:), pointer :: sfcMassBalTimeSeries, sfcAirTempTimeSeries, betaTimeSeries, basalHeatFluxTimeSeries, marineBasalMassBalTimeSeries
+ real (kind=RKIND), dimension(:), pointer :: sfcMassBal, sfcAirTemp, beta, basalHeatFlux, marineBasalMassBal
integer :: iYear, ierr, nCells
@@ -105,11 +105,13 @@
sfcAirTempTimeSeries => grid % sfcAirTempTimeSeries % array
betaTimeSeries => grid % betaTimeSeries % array
basalHeatFluxTimeSeries => grid % basalHeatFluxTimeSeries % array
+ marineBasalMassBalTimeSeries => grid % marineBasalMassBalTimeSeries % array
sfcMassBal => grid % sfcMassBal % array
sfcAirTemp => grid % sfcAirTemp % array
beta => grid % beta % array
basalHeatFlux => grid % basalHeatFlux % array
+ marineBasalMassBal => grid % marineBasalMassBal % array
nCells = grid % nCells
@@ -127,6 +129,7 @@
call assign_annual_forcing_to_1D_array(sfcAirTemp, sfcAirTempTimeSeries, iYear, nCells)
call assign_annual_forcing_to_1D_array(beta, betaTimeSeries, iYear, nCells)
call assign_annual_forcing_to_1D_array(basalHeatFlux, basalHeatFluxTimeSeries, iYear, nCells)
+ call assign_annual_forcing_to_1D_array(marineBasalMassBal, marineBasalMassBalTimeSeries, iYear, nCells)
end subroutine land_ice_assign_annual_forcing!}}}
@@ -197,7 +200,6 @@
end subroutine assign_annual_forcing_to_1D_array
-
!***********************************************************************
end module land_ice_annual_forcing
Modified: branches/land_ice_projects/implement_core/src/core_land_ice/mpas_land_ice_lifev.F
===================================================================
--- branches/land_ice_projects/implement_core/src/core_land_ice/mpas_land_ice_lifev.F        2013-01-30 22:34:57 UTC (rev 2399)
+++ branches/land_ice_projects/implement_core/src/core_land_ice/mpas_land_ice_lifev.F        2013-01-31 03:40:57 UTC (rev 2400)
@@ -104,7 +104,8 @@
!call lifeV and set the grid of the velocity solver
call velocity_solver_init_mpi(domain % dminfo % comm)
#else
- write(*,*) "LifeV needed to run L1L2, FO and Stokes models"
+ write(*,*) "LifeV needed to run L1L2, FO and Stokes models" ! stdout
+ write(0,*) "LifeV needed to run L1L2, FO and Stokes models" ! stderr
err = 1
return
#endif
@@ -163,6 +164,7 @@
call phg_init(domain % dminfo % comm)
#else
write(*,*) "PHG needed to run Stokes model"
+ write(0,*) "PHG needed to run Stokes model"
err = 1
return
#endif
@@ -341,6 +343,7 @@
call mpas_timer_stop("velocity_solver_set_grid_data")
#else
write(*,*) "LifeV needed to run L1L2, FO and Stokes models"
+ write(0,*) "LifeV needed to run L1L2, FO and Stokes models"
err = 1
return
#endif
Modified: branches/land_ice_projects/implement_core/src/core_land_ice/mpas_land_ice_tendency.F
===================================================================
--- branches/land_ice_projects/implement_core/src/core_land_ice/mpas_land_ice_tendency.F        2013-01-30 22:34:57 UTC (rev 2399)
+++ branches/land_ice_projects/implement_core/src/core_land_ice/mpas_land_ice_tendency.F        2013-01-31 03:40:57 UTC (rev 2400)
@@ -112,7 +112,10 @@
! local variables
!
!-----------------------------------------------------------------
+ integer nVertLevels
+ nVertLevels = mesh % nVertLevels
+
! 0 tendency
layerThickness_tend = 0.0
@@ -156,11 +159,23 @@
! Do nothing - don't add the MB
case default
! Add surface mass balance to tendency
- ! NOTE: Need to decide where to put this for layer by layer thickness tendency, and how to deal with negative SMB that eliminates top layer or all ice (check for negative thickness?)
- layerThickness_tend(1,:) = layerThickness_tend(1,:) + mesh % sfcMassBal % array ! (tendency in meters per year)
- ! THIS MIGHT RESULT IN NEGATIVE LAYER THICKNESS!
+ ! TODO: Need to decide where to put this for layer by layer thickness tendency, and how to deal with negative SMB that eliminates top layer or all ice (check for negative thickness?)
+ layerThickness_tend(1,:) = layerThickness_tend(1,:) + mesh % sfcMassBal % array ! (tendency in meters per year)
+ ! TODO THIS MIGHT RESULT IN NEGATIVE LAYER THICKNESS!
- ! Do basal mass balance too?
+ ! Add basal mass balance to tendency
+ ! Apply marineBasalMassBal to floating ice only. Apply groundedBasalMassBal to grounded ice only.
+ ! TODO: more complicated treatment at GL?
+ ! Assign 0.0 values to the marineBasalMassBal where the ice is grounded.
+ ! It's ok to overwrite the values with 0's here, because each time step
+ ! we get a fresh copy of the array from the annual_forcing subroutine.
+ where ( MASK_IS_GROUNDED(state % cellMask % array) )
+ mesh % marineBasalMassBal % array = 0.0
+ end where
+ layerThickness_tend(nVertLevels,:) = layerThickness_tend(nVertLevels,:) &
+ + mesh % marineBasalMassBal % array ! (tendency in meters per year)
+ ! TODO Add in grounded ice basal mass balance once temperature diffusion is calculated
+ ! TODO THIS MIGHT RESULT IN NEGATIVE LAYER THICKNESS!
end select
@@ -346,7 +361,7 @@
!-----------------------------------------------------------------
integer :: iCell, iLevel, nCells, nVertLevels
real (kind=RKIND), dimension(:), pointer :: thickness, upperSurface, &
- lowerSurface, bedTopography, layerThicknessFractions
+ lowerSurface, bedTopography, layerThicknessFractions, beta
integer, dimension(:), pointer :: cellMask
real (kind=RKIND), dimension(:,:), pointer :: layerThickness
@@ -359,15 +374,23 @@
lowerSurface => state % lowerSurface % array
bedTopography => state % bedTopography % array
layerThickness => state % layerThickness % array
+ beta => mesh % beta % array
cellMask => state % cellMask % array
! Calculate masks - needs to happen before calculating lower surface so we know where the ice is floating
call land_ice_calculate_mask(mesh, state, err)
- ! no ice shelves -- lower surface same as bed topography
- ! Eventually this should consider floatation!
+ ! Update beta before the velocity solve occurs, now that we have the new state and its mask.
+ ! It's ok to overwrite the beta values with 0's here, because each time step
+ ! we get a fresh copy of the array from the annual_forcing subroutine.
+ ! Note: some velocity solvers may do this on their own, but we are doing it here for completeness.
where ( MASK_IS_FLOATING(cellMask) )
+ beta = 0.0
+ end where
+
+ ! Lower surface is based on floatation for floating ice. For grounded ice it is the bed.
+ where ( MASK_IS_FLOATING(cellMask) )
lowerSurface = config_sealevel - thickness * (config_ice_density / config_ocean_density)
elsewhere
lowerSurface = bedTopography
</font>
</pre>