MOM6
MOM_ice_shelf_initialize.F90
Go to the documentation of this file.
1 !> Initialize ice shelf variables
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
6 use mom_grid, only : ocean_grid_type
8 use mom_io, only: mom_read_data, file_exists, slasher
9 use mom_error_handler, only : mom_error, mom_mesg, fatal, warning, is_root_pe
12 
13 implicit none ; private
14 
15 #include <MOM_memory.h>
16 
17 !MJHpublic initialize_ice_shelf_boundary, initialize_ice_thickness
19 
20 ! A note on unit descriptions in comments: MOM6 uses units that can be rescaled for dimensional
21 ! consistency testing. These are noted in comments with units like Z, H, L, and T, along with
22 ! their mks counterparts with notation like "a velocity [Z T-1 ~> m s-1]". If the units
23 ! vary with the Boussinesq approximation, the Boussinesq variant is given first.
24 
25 contains
26 
27 !> Initialize ice shelf thickness
28 subroutine initialize_ice_thickness(h_shelf, area_shelf_h, hmask, G, US, PF)
29  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
30  real, dimension(SZDI_(G),SZDJ_(G)), &
31  intent(inout) :: h_shelf !< The ice shelf thickness [Z ~> m].
32  real, dimension(SZDI_(G),SZDJ_(G)), &
33  intent(inout) :: area_shelf_h !< The area per cell covered by the ice shelf [m2].
34  real, dimension(SZDI_(G),SZDJ_(G)), &
35  intent(inout) :: hmask !< A mask indicating which tracer points are
36  !! partly or fully covered by an ice-shelf
37  type(unit_scale_type), intent(in) :: us !< A structure containing unit conversion factors
38  type(param_file_type), intent(in) :: pf !< A structure to parse for run-time parameters
39 
40  integer :: i, j
41  character(len=40) :: mdl = "initialize_ice_thickness" ! This subroutine's name.
42  character(len=200) :: config
43 
44  call get_param(pf, mdl, "ICE_PROFILE_CONFIG", config, &
45  "This specifies how the initial ice profile is specified. "//&
46  "Valid values are: CHANNEL, FILE, and USER.", &
47  fail_if_missing=.true.)
48 
49  select case ( trim(config) )
50  case ("CHANNEL"); call initialize_ice_thickness_channel (h_shelf, area_shelf_h, hmask, g, us, pf)
51  case ("FILE"); call initialize_ice_thickness_from_file (h_shelf, area_shelf_h, hmask, g, us, pf)
52  case ("USER"); call user_init_ice_thickness (h_shelf, area_shelf_h, hmask, g, us, pf)
53  case default ; call mom_error(fatal,"MOM_initialize: "// &
54  "Unrecognized ice profile setup "//trim(config))
55  end select
56 
57 end subroutine initialize_ice_thickness
58 
59 !> Initialize ice shelf thickness from file
60 subroutine initialize_ice_thickness_from_file(h_shelf, area_shelf_h, hmask, G, US, PF)
61  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
62  real, dimension(SZDI_(G),SZDJ_(G)), &
63  intent(inout) :: h_shelf !< The ice shelf thickness [m].
64  real, dimension(SZDI_(G),SZDJ_(G)), &
65  intent(inout) :: area_shelf_h !< The area per cell covered by the ice shelf [m2].
66  real, dimension(SZDI_(G),SZDJ_(G)), &
67  intent(inout) :: hmask !< A mask indicating which tracer points are
68  !! partly or fully covered by an ice-shelf
69  type(unit_scale_type), intent(in) :: US !< A structure containing unit conversion factors
70  type(param_file_type), intent(in) :: PF !< A structure to parse for run-time parameters
71 
72  ! This subroutine reads ice thickness and area from a file and puts it into
73  ! h_shelf and area_shelf_h in m (and dimensionless) and updates hmask
74  character(len=200) :: filename,thickness_file,inputdir ! Strings for file/path
75  character(len=200) :: thickness_varname, area_varname ! Variable name in file
76  character(len=40) :: mdl = "initialize_ice_thickness_from_file" ! This subroutine's name.
77  integer :: i, j, isc, jsc, iec, jec
78  real :: len_sidestress, mask, udh
79 
80  call mom_mesg(" MOM_ice_shelf_init_profile.F90, initialize_thickness_from_file: reading thickness")
81 
82  call get_param(pf, mdl, "INPUTDIR", inputdir, default=".")
83  inputdir = slasher(inputdir)
84  call get_param(pf, mdl, "ICE_THICKNESS_FILE", thickness_file, &
85  "The file from which the bathymetry is read.", &
86  default="ice_shelf_h.nc")
87  call get_param(pf, mdl, "LEN_SIDE_STRESS", len_sidestress, &
88  "position past which shelf sides are stress free.", &
89  default=0.0, units="axis_units")
90 
91  filename = trim(inputdir)//trim(thickness_file)
92  call log_param(pf, mdl, "INPUTDIR/THICKNESS_FILE", filename)
93  call get_param(pf, mdl, "ICE_THICKNESS_VARNAME", thickness_varname, &
94  "The name of the thickness variable in ICE_THICKNESS_FILE.", &
95  default="h_shelf")
96  call get_param(pf, mdl, "ICE_AREA_VARNAME", area_varname, &
97  "The name of the area variable in ICE_THICKNESS_FILE.", &
98  default="area_shelf_h")
99 
100  if (.not.file_exists(filename, g%Domain)) call mom_error(fatal, &
101  " initialize_topography_from_file: Unable to open "//trim(filename))
102 
103  call mom_read_data(filename, trim(thickness_varname), h_shelf, g%Domain, scale=us%m_to_Z)
104  call mom_read_data(filename,trim(area_varname),area_shelf_h,g%Domain)
105 
106 ! call get_param(PF, mdl, "ICE_BOUNDARY_CONFIG", config, &
107 ! "This specifies how the ice domain boundary is specified", &
108 ! fail_if_missing=.true.)
109 
110  isc = g%isc ; jsc = g%jsc ; iec = g%iec ; jec = g%jec
111 
112  do j=jsc,jec
113  do i=isc,iec
114 
115  ! taper ice shelf in area where there is no sidestress -
116  ! but do not interfere with hmask
117 
118  if ((g%geoLonCv(i,j) > len_sidestress).and. &
119  (len_sidestress > 0.)) then
120  udh = exp(-(g%geoLonCv(i,j)-len_sidestress)/5.0) * h_shelf(i,j)
121  if (udh <= 25.0) then
122  h_shelf(i,j) = 0.0
123  area_shelf_h(i,j) = 0.0
124  else
125  h_shelf(i,j) = udh
126  endif
127  endif
128 
129  ! update thickness mask
130 
131  if (area_shelf_h(i,j) >= us%L_to_m**2*g%areaT(i,j)) then
132  hmask(i,j) = 1.
133  elseif (area_shelf_h(i,j) == 0.0) then
134  hmask(i,j) = 0.
135  elseif ((area_shelf_h(i,j) > 0) .and. (area_shelf_h(i,j) <= us%L_to_m**2*g%areaT(i,j))) then
136  hmask(i,j) = 2.
137  else
138  call mom_error(fatal,mdl// " AREA IN CELL OUT OF RANGE")
139  endif
140  enddo
141  enddo
142 
143 
145 
146 !> Initialize ice shelf thickness for a channel configuration
147 subroutine initialize_ice_thickness_channel(h_shelf, area_shelf_h, hmask, G, US, PF)
148  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
149  real, dimension(SZDI_(G),SZDJ_(G)), &
150  intent(inout) :: h_shelf !< The ice shelf thickness [Z ~> m].
151  real, dimension(SZDI_(G),SZDJ_(G)), &
152  intent(inout) :: area_shelf_h !< The area per cell covered by the ice shelf [m2].
153  real, dimension(SZDI_(G),SZDJ_(G)), &
154  intent(inout) :: hmask !< A mask indicating which tracer points are
155  !! partly or fully covered by an ice-shelf
156  type(unit_scale_type), intent(in) :: US !< A structure containing unit conversion factors
157  type(param_file_type), intent(in) :: PF !< A structure to parse for run-time parameters
158 
159  character(len=40) :: mdl = "initialize_ice_shelf_thickness_channel" ! This subroutine's name.
160  real :: max_draft, min_draft, flat_shelf_width, c1, slope_pos
161  real :: edge_pos, shelf_slope_scale, Rho_ocean
162  integer :: i, j, jsc, jec, jsd, jed, jedg, nyh, isc, iec, isd, ied
163  integer :: j_off
164 
165  jsc = g%jsc ; jec = g%jec ; isc = g%isc ; iec = g%iec
166  jsd = g%jsd ; jed = g%jed ; isd = g%isd ; ied = g%ied
167  nyh = g%domain%njhalo ; jedg = g%domain%njglobal+nyh
168  j_off = g%jdg_offset
169 
170  call mom_mesg(mdl//": setting thickness")
171 
172  call get_param(pf, mdl, "SHELF_MAX_DRAFT", max_draft, &
173  units="m", default=1.0, scale=us%m_to_Z)
174  call get_param(pf, mdl, "SHELF_MIN_DRAFT", min_draft, &
175  units="m", default=1.0, scale=us%m_to_Z)
176  call get_param(pf, mdl, "FLAT_SHELF_WIDTH", flat_shelf_width, &
177  units="axis_units", default=0.0)
178  call get_param(pf, mdl, "SHELF_SLOPE_SCALE", shelf_slope_scale, &
179  units="axis_units", default=0.0)
180  call get_param(pf, mdl, "SHELF_EDGE_POS_0", edge_pos, &
181  units="axis_units", default=0.0)
182 ! call get_param(param_file, mdl, "RHO_0", Rho_ocean, &
183 ! "The mean ocean density used with BOUSSINESQ true to "//&
184 ! "calculate accelerations and the mass for conservation "//&
185 ! "properties, or with BOUSSINSEQ false to convert some "//&
186 ! "parameters from vertical units of m to kg m-2.", &
187 ! units="kg m-3", default=1035.0, scale=US%Z_to_m)
188 
189  slope_pos = edge_pos - flat_shelf_width
190  c1 = 0.0 ; if (shelf_slope_scale > 0.0) c1 = 1.0 / shelf_slope_scale
191 
192 
193  do j=g%jsd,g%jed
194 
195  if (((j+j_off) <= jedg) .AND. ((j+j_off) >= nyh+1)) then
196 
197  do i=g%isc,g%iec
198 
199  if ((j >= jsc) .and. (j <= jec)) then
200 
201  if (g%geoLonCu(i-1,j) >= edge_pos) then
202  ! Everything past the edge is open ocean.
203 ! mass_shelf(i,j) = 0.0
204  area_shelf_h(i,j) = 0.0
205  hmask(i,j) = 0.0
206  h_shelf(i,j) = 0.0
207  else
208  if (g%geoLonCu(i,j) > edge_pos) then
209  area_shelf_h(i,j) = us%L_to_m**2*g%areaT(i,j) * (edge_pos - g%geoLonCu(i-1,j)) / &
210  (g%geoLonCu(i,j) - g%geoLonCu(i-1,j))
211  hmask(i,j) = 2.0
212  else
213  area_shelf_h(i,j) = us%L_to_m**2*g%areaT(i,j)
214  hmask(i,j) = 1.0
215  endif
216 
217  if (g%geoLonT(i,j) > slope_pos) then
218  h_shelf(i,j) = min_draft
219 ! mass_shelf(i,j) = Rho_ocean * min_draft
220  else
221 ! mass_shelf(i,j) = Rho_ocean * (min_draft + &
222 ! (CS%max_draft - CS%min_draft) * &
223 ! min(1.0, (c1*(slope_pos - G%geoLonT(i,j)))**2) )
224  h_shelf(i,j) = (min_draft + &
225  (max_draft - min_draft) * &
226  min(1.0, (c1*(slope_pos - g%geoLonT(i,j)))**2) )
227  endif
228 
229  endif
230  endif
231 
232  if ((i+g%idg_offset) == g%domain%nihalo+1) then
233  hmask(i-1,j) = 3.0
234  endif
235 
236  enddo
237  endif ; enddo
238 
240 
241 !BEGIN MJH
242 ! subroutine initialize_ice_shelf_boundary(u_face_mask_bdry, v_face_mask_bdry, &
243 ! u_flux_bdry_val, v_flux_bdry_val, u_bdry_val, v_bdry_val, h_bdry_val, &
244 ! hmask, G, PF )
245 
246 ! type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
247 ! real, dimension(SZIB_(G),SZJ_(G)), &
248 ! intent(inout) :: u_face_mask_bdry !< A boundary-type mask at C-grid u faces
249 ! real, dimension(SZIB_(G),SZJ_(G)), &
250 ! intent(inout) :: u_flux_bdry_val !< The boundary thickness flux through
251  !! C-grid u faces [m2 s-1].
252 ! real, dimension(SZI_(G),SZJB_(G)), &
253 ! intent(inout) :: v_face_mask_bdry !< A boundary-type mask at C-grid v faces
254 ! real, dimension(SZI_(G),SZJB_(G)), &
255 ! intent(inout) :: v_flux_bdry_val !< The boundary thickness flux through
256  !! C-grid v faces [m2 s-1].
257 ! real, dimension(SZIB_(G),SZJB_(G)), &
258 ! intent(inout) :: u_bdry_val !< The zonal ice shelf velocity at open
259  !! boundary vertices [m yr-1].
260 ! real, dimension(SZIB_(G),SZJB_(G)), &
261 ! intent(inout) :: v_bdry_val !< The meridional ice shelf velocity at open
262  !! boundary vertices [m yr-1].
263 ! real, dimension(SZDI_(G),SZDJ_(G)), &
264 ! intent(inout) :: h_bdry_val !< The ice shelf thickness at open boundaries
265 ! real, dimension(SZDI_(G),SZDJ_(G)), &
266 ! intent(inout) :: hmask !< A mask indicating which tracer points are
267 ! !! partly or fully covered by an ice-shelf
268 ! type(param_file_type), intent(in) :: PF !< A structure to parse for run-time parameters
269 
270 ! character(len=40) :: mdl = "initialize_ice_shelf_boundary" ! This subroutine's name.
271 ! character(len=200) :: config
272 ! logical flux_bdry
273 
274 ! call get_param(PF, mdl, "ICE_BOUNDARY_CONFIG", config, &
275 ! "This specifies how the ice domain boundary is specified. "//&
276 ! "valid values include CHANNEL, FILE and USER.", &
277 ! fail_if_missing=.true.)
278 ! call get_param(PF, mdl, "ICE_BOUNDARY_FLUX_CONDITION", flux_bdry, &
279 ! "This specifies whether mass input is a dirichlet or "//&
280 ! "flux condition", default=.true.)
281 
282 ! select case ( trim(config) )
283 ! case ("CHANNEL")
284 ! call initialize_ice_shelf_boundary_channel(u_face_mask_bdry, &
285 ! v_face_mask_bdry, u_flux_bdry_val, v_flux_bdry_val, &
286 ! u_bdry_val, v_bdry_val, h_bdry_val, hmask, G, &
287 ! flux_bdry, PF)
288 ! case ("FILE"); call MOM_error(FATAL,"MOM_initialize: "// &
289 ! "Unrecognized topography setup "//trim(config))
290 ! case ("USER"); call MOM_error(FATAL,"MOM_initialize: "// &
291 ! "Unrecognized topography setup "//trim(config))
292 ! case default ; call MOM_error(FATAL,"MOM_initialize: "// &
293 ! "Unrecognized topography setup "//trim(config))
294 ! end select
295 
296 ! end subroutine initialize_ice_shelf_boundary
297 
298 ! subroutine initialize_ice_shelf_boundary_channel(u_face_mask_bdry, v_face_mask_bdry, &
299 ! u_flux_bdry_val, v_flux_bdry_val, u_bdry_val, v_bdry_val, h_bdry_val, &
300 ! hmask, G, flux_bdry, PF )
301 
302 ! type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
303 ! real, dimension(SZIB_(G),SZJ_(G)), &
304 ! intent(inout) :: u_face_mask_bdry !< A boundary-type mask at C-grid u faces
305 ! real, dimension(SZIB_(G),SZJ_(G)), &
306 ! intent(inout) :: u_flux_bdry_val !< The boundary thickness flux through
307  !! C-grid u faces [m2 s-1].
308 ! real, dimension(SZI_(G),SZJB_(G)), &
309 ! intent(inout) :: v_face_mask_bdry !< A boundary-type mask at C-grid v faces
310 ! real, dimension(SZI_(G),SZJB_(G)), &
311 ! intent(inout) :: v_flux_bdry_val !< The boundary thickness flux through
312  !! C-grid v faces [m2 s-1].
313 ! real, dimension(SZIB_(G),SZJB_(G)), &
314 ! intent(inout) :: u_bdry_val !< The zonal ice shelf velocity at open
315  !! boundary vertices [m yr-1].
316 ! real, dimension(SZIB_(G),SZJB_(G)), &
317 ! intent(inout) :: v_bdry_val !< The meridional ice shelf velocity at open
318  !! boundary vertices [m yr-1].
319 ! real, dimension(SZDI_(G),SZDJ_(G)), &
320 ! intent(inout) :: h_bdry_val !< The ice shelf thickness at open boundaries
321 ! real, dimension(SZDI_(G),SZDJ_(G)), &
322 ! intent(inout) :: hmask !< A mask indicating which tracer points are
323 ! !! partly or fully covered by an ice-shelf
324 ! logical, intent(in) :: flux_bdry !< If true, use mass fluxes as the boundary value.
325 ! type(param_file_type), intent(in) :: PF !< A structure to parse for run-time parameters
326 
327 ! character(len=40) :: mdl = "initialize_ice_shelf_boundary_channel" ! This subroutine's name.
328 ! integer :: i, j, isd, jsd, is, js, iegq, jegq, giec, gjec, gisc, gjsc, isc, jsc, iec, jec, ied, jed
329 ! real :: lenlat, input_thick, input_flux, len_stress
330 
331 ! call get_param(PF, mdl, "LENLAT", lenlat, fail_if_missing=.true.)
332 
333 ! call get_param(PF, mdl, "INPUT_FLUX_ICE_SHELF", input_flux, &
334 ! "volume flux at upstream boundary", &
335 ! units="m2 s-1", default=0.)
336 ! call get_param(PF, mdl, "INPUT_THICK_ICE_SHELF", input_thick, &
337 ! "flux thickness at upstream boundary", &
338 ! units="m", default=1000.)
339 ! call get_param(PF, mdl, "LEN_SIDE_STRESS", len_stress, &
340 ! "maximum position of no-flow condition in along-flow direction", &
341 ! units="km", default=0.)
342 
343 ! call MOM_mesg(mdl//": setting boundary")
344 
345 ! isd = G%isd ; ied = G%ied
346 ! jsd = G%jsd ; jed = G%jed
347 ! isc = G%isc ; jsc = G%jsc ; iec = G%iec ; jec = G%jec
348 ! gisc = G%Domain%nihalo ; gjsc = G%Domain%njhalo
349 ! giec = G%Domain%niglobal+gisc ; gjec = G%Domain%njglobal+gjsc
350 
351 ! do j=jsd,jed
352 ! do i=isd,ied
353 
354 ! ! upstream boundary - set either dirichlet or flux condition
355 
356 ! if ((i+G%idg_offset) == G%domain%nihalo+1) then
357 ! if (flux_bdry) then
358 ! u_face_mask_bdry(i-1,j) = 4.0
359 ! u_flux_bdry_val(i-1,j) = input_flux
360 ! else
361 ! hmask(i-1,j) = 3.0
362 ! h_bdry_val(i-1,j) = input_thick
363 ! u_face_mask_bdry(i-1,j) = 3.0
364 ! u_bdry_val(i-1,j-1) = (1 - ((G%geoLatBu(i-1,j-1) - 0.5*lenlat)*2./lenlat)**2) * &
365 ! 1.5 * input_flux / input_thick
366 ! u_bdry_val(i-1,j) = (1 - ((G%geoLatBu(i-1,j) - 0.5*lenlat)*2./lenlat)**2) * &
367 ! 1.5 * input_flux / input_thick
368 ! endif
369 ! endif
370 
371 ! ! side boundaries: no flow
372 
373 ! if (G%jdg_offset+j == gjsc+1) then !bot boundary
374 ! if (len_stress == 0. .OR. G%geoLonCv(i,j-1) <= len_stress) then
375 ! v_face_mask_bdry(i,j-1) = 0.
376 ! else
377 ! v_face_mask_bdry(i,j-1) = 1.
378 ! endif
379 ! elseif (G%jdg_offset+j == gjec) then !top boundary
380 ! if (len_stress == 0. .OR. G%geoLonCv(i,j-1) <= len_stress) then
381 ! v_face_mask_bdry(i,j) = 0.
382 ! else
383 ! v_face_mask_bdry(i,j) = 1.
384 ! endif
385 ! endif
386 
387 ! ! downstream boundary - CFBC
388 
389 ! if (i+G%idg_offset == giec) then
390 ! u_face_mask_bdry(i,j) = 2.0
391 ! endif
392 
393 ! enddo
394 ! enddo
395 
396 !END MJH end subroutine initialize_ice_shelf_boundary_channel
397 
398 end module mom_ice_shelf_initialize
user_shelf_init::user_init_ice_thickness
subroutine, public user_init_ice_thickness(h_shelf, area_shelf_h, hmask, G, US, param_file)
This subroutine updates the ice shelf thickness, as specified by user-provided code.
Definition: user_shelf_init.F90:104
mom_ice_shelf_initialize::initialize_ice_thickness_from_file
subroutine initialize_ice_thickness_from_file(h_shelf, area_shelf_h, hmask, G, US, PF)
Initialize ice shelf thickness from file.
Definition: MOM_ice_shelf_initialize.F90:61
mom_ice_shelf_initialize
Initialize ice shelf variables.
Definition: MOM_ice_shelf_initialize.F90:2
mom_error_handler::mom_mesg
subroutine, public mom_mesg(message, verb, all_print)
This provides a convenient interface for writing an informative comment.
Definition: MOM_error_handler.F90:53
mom_file_parser::param_file_type
A structure that can be parsed to read and document run-time parameters.
Definition: MOM_file_parser.F90:54
mom_file_parser::get_param
An overloaded interface to read and log the values of various types of parameters.
Definition: MOM_file_parser.F90:102
mom_io
This module contains I/O framework code.
Definition: MOM_io.F90:2
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
mom_ice_shelf_initialize::initialize_ice_thickness_channel
subroutine initialize_ice_thickness_channel(h_shelf, area_shelf_h, hmask, G, US, PF)
Initialize ice shelf thickness for a channel configuration.
Definition: MOM_ice_shelf_initialize.F90:148
mom_ice_shelf_initialize::initialize_ice_thickness
subroutine, public initialize_ice_thickness(h_shelf, area_shelf_h, hmask, G, US, PF)
Initialize ice shelf thickness.
Definition: MOM_ice_shelf_initialize.F90:29
mom_io::mom_read_data
Read a data field from a file.
Definition: MOM_io.F90:74
mom_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_unit_scaling
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
Definition: MOM_unit_scaling.F90:2
mom_error_handler::is_root_pe
logical function, public is_root_pe()
This returns .true. if the current PE is the root PE.
Definition: MOM_error_handler.F90:44
mom_io::file_exists
Indicate whether a file exists, perhaps with domain decomposition.
Definition: MOM_io.F90:68
mom_error_handler::mom_error
subroutine, public mom_error(level, message, all_print)
This provides a convenient interface for writing an mpp_error message with run-time filter based on a...
Definition: MOM_error_handler.F90:72
user_shelf_init
This module specifies the initial values and evolving properties of the MOM6 ice shelf,...
Definition: user_shelf_init.F90:3
mom_file_parser::log_param
An overloaded interface to log the values of various types of parameters.
Definition: MOM_file_parser.F90:96
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26
mom_file_parser::read_param
An overloaded interface to read various types of parameters.
Definition: MOM_file_parser.F90:90