This appendix details the subroutine calls provided by GLINT, and their arguments. Note that where a type is given as real(rk), this indicates that the kind of the real type is specified by the value of the parameter rk, which may be altered at compile-time (see appropriate other documentation for details).
Purpose To initialise the ice model, and load in all relevant parameter files.
subroutine initialise_glint(params,lats,longs,paramfile,time_step)
|
| Mandatory
| |||
| params | type(glint_params) | intent(inout) | Ice model to be configured |
| lats(:) | real(rk) | intent(in) | latitudinal location of grid-points in global data (given in ∘N) |
| longs(:) | real(rk) | intent(in) | longitudinal location of grid-points in global data (given in ∘E) |
| paramfile | character(*) | intent(in) | name of top-level parameter file |
| time_step | integer | intent(in) | Intended calling time-step (hours) |
| Optional
| |||
| latb(:) | real(rk) | intent(in) | Latiudinal locations of grid-box boundaries (degrees). This array has one more element than lats. |
| lonb(:) | real(rk) | intent(in) | Longitudinal locations of grid-box boundaries (degrees). This array has one more element than longs. |
| orog(:,:) | real(rk) | intent(out) | The initial orography (m). |
| albedo | real(rk) | intent(out) | The initial ice albedo field |
| ice_frac | real(rk) | intent(out) | The initial ice fraction |
| orog_lats | real(rk) | intent(in) | Latitudinal location of gridpoints for global orography output |
| orog_longs | real(rk) | intent(in) | Longitudinal location of gridpoints for global orography output |
| orog_latb | real(rk) | intent(in) | Locations of the latitudinal boundaries of the grid-boxes (orography) |
| orog_lonb | real(rk) | intent(in) | Locations of the longitudinal boundaries of the grid-boxes (orography) |
| output_flag | logical | intent(out) | Set to show outputs have been updated (provided for consistency with main glint subroutine). |
Purpose To perform temporal averaging of input fields, and, if necessary, down-scale those fields onto local projections and perform an ice model time-step. Output files may be appended to, and if optional arguments used, fields made available for feedback.
subroutine glint(params,time,temp,precip,orog)
|
| Mandatory
| |||
| params | type(glint_params) | intent(inout) | parameters for this run |
| time | integer | intent(in) | Current model time (hours) |
| temp(:,:) | real(rk) | intent(in) | Surface temperature field (∘C) |
| precip(:,:) | real(rk) | intent(in) | Precipitation field (mm/s) |
| orog(:,:) | real(rk) | intent(in) | Global orography (m) |
| Optional
| |||
| zonwind(:,:) | real(rk) | intent(in) | Zonal component of the wind field (ms-1) |
| merwind(:,:) | real(rk) | intent(in) | Meridional component of the wind field (ms-1) |
| output_flag | logical | intent(out) | Set to show new output fields have been calculated after an ice-model time-step. If this flag is not set, the output fields retain their values at input. |
| orog_out(:,:) | real(rk) | intent(inout) | Output orography (m) |
| albedo(:,:) | real(rk) | intent(inout) | Surface albedo |
| ice_frac(:,:) | real(rk) | intent(inout) | Fractional ice coverage |
| water_in(:,:) | real(rk) | intent(inout) | The input fresh-water flux (mm, over ice time-step). Essentially precip, but provided for consistency. |
| water_out(:,:) | real(rk) | intent(inout) | The output fresh-water flux (mm, over ice time-step). This is simply the ablation calculated by the model, scaled up to the global grid. It is up to the global model to then deal with it (route it to the oceans, land scheme, etc.) Note that the precipitation fed to the model but which doesn’t get incorporated into the ice sheet because it falls over the sea is returned in this field. |
| total_water_in | real(rk) | intent(inout) | Area-integrated water flux in (kg) |
| total_water_out | real(rk) | intent(inout) | Area-integrated water flux out (kg) |
| ice_volume | real(rk) | intent(inout) | Total ice volume (m3) |
Example interpretation of output fields Consider a particular point, (i,j) in the global domain. Suppose value returned by glint_coverage_map for this point is 0.7, and the output fields have these values:
orog_out(i,j) = 200.0
albedo(i,j) = 0.4 ice_frac(i,j) = 0.5 |
What does this mean? Well, the ice model covers 70% of the grid-box, and in that part the mean surface elevation is 200 m. Of the part covered by the ice model, half is actually covered by ice. Thus, 35% (0.5 × 0.7) of the global grid-box is covered by ice, and the ice has an mean albedo of 40%. The model makes no suggestion for the albedo or elevation of the other 65% of the grid-box. Currently, ice albedo is a constant that may be changed in the appropriate configuration file, but this output field is provided against the possibility that the model may be extended at some point to include a model of ice albedo.
Purpose To perform general tidying-up operations, close files, etc.
subroutine end_glint(params)
|
| params | type(glint_params) | intent(inout) | Ice model paramters |
Purpose To obtain a map of fractional coverage of global grid-boxes by the ice model instances. The function returns a value indicating success, or giving error information.
Type, name and mandatory arguments
integer function glint_coverage_map(params,coverage,cov_orog)
|
| params | type(glint_params) | intent(in) | Ice model parameters |
| coverage(:,:) | real(rk) | intent(out) | Coverage map (all fields except orography) |
| cov_orog(:,:) | real(rk) | intent(out) | Coverage map (orography) |
| Value | Meaning |
| 0 | Coverage maps have been returned successfully |
| 1 | Coverage maps not yet calculated; must call initialise_glint first |
| 2 | Arrays coverage or cov_orog are the wrong size |