!||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| module forcing_ws !BOP ! !MODULE: forcing_ws ! ! !DESCRIPTION: ! Contains routines and variables used for determining the ! surface wind stress. ! ! !REVISION HISTORY: ! SVN:$Id: forcing_ws.F90 12674 2008-10-31 22:21:32Z njn01 $ ! ! !USES: use kinds_mod use constants use blocks use distribution use domain use io use forcing_tools use time_management use grid use exit_mod implicit none private save ! !PUBLIC MEMBER FUNCTIONS: public :: init_ws, & set_ws ! !PUBLIC DATA MEMBERS: real (r8), public :: &! needed by restart ws_interp_last ! time when last interpolation was done !EOP !BOC !----------------------------------------------------------------------- ! ! module variables ! !----------------------------------------------------------------------- real (r8), allocatable, dimension(:,:,:,:,:) :: & SMF_DATA ! surface momentum flux data at multiple times that ! may be interpolated in time to get SMF real (r8), dimension(12) :: & ws_data_time ! time (hours) corresponding to surface fluxes ! in SMF_DATA real (r8), dimension(20) :: & ws_data_renorm ! factors for converting to model units real (r8) :: & ws_data_inc, &! time increment between values of forcing data ws_data_next, &! time for next value of forcing data needed ws_data_update, &! time new forcing value needs to be added to interpolation set ws_interp_inc, &! time increment between interpolation ws_interp_next ! time when next interpolation will be done integer (int_kind) :: & ws_interp_order, &! order of temporal interpolation ws_data_time_min_loc ! index of first time index of SMF_DATA ! to use for interpolating forcing character (char_len) :: & ws_filename, &! name of file conainting forcing data ws_file_fmt, &! format (bin or nc) for forcing file ws_interp_freq, &! keyword for period of temporal interpolation ws_interp_type, &! ws_data_label, &! name of data to be read ws_formulation ! formulation to use to compute wind stress character (char_len), public :: & ws_data_type ! keyword for period of the forcing data character (char_len), dimension(:), allocatable :: & ws_data_names ! names of each data field integer (int_kind), dimension(:), allocatable :: & ws_bndy_loc, &! location and field type for ghost cell ws_bndy_type ! updates !EOC !*********************************************************************** contains !*********************************************************************** !BOP ! !IROUTINE: init_ws ! !INTERFACE: subroutine init_ws(SMF,SMFT,lsmft_avail) ! !DESCRIPTION: ! Initializes wind stress forcing by either calculating or reading ! in the wind stress. Also performs initial book-keeping concerning ! when new data is needed for the temporal interpolation and ! when forcing will need to be updated. ! ! !REVISION HISTORY: ! same as module ! !OUTPUT PARAMETERS: logical (log_kind), intent(out) :: & lsmft_avail ! flag is true if momentum fluxes available ! at T points (to prevent later averaging) real (r8), dimension(nx_block,ny_block,2,max_blocks_clinic), & intent(out) :: & SMF, &! surface momentum fluxes at current timestep SMFT ! surface momentum fluxes at T points !EOP !BOC !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- integer (int_kind) :: & n, &! loop index iblock, &! block index nml_error ! namelist i/o error flag character (char_len) :: & forcing_filename ! full name of forcing input file namelist /forcing_ws_nml/ ws_data_type, ws_data_inc, & ws_interp_type, ws_interp_freq, & ws_interp_inc, ws_filename, & ws_file_fmt, ws_data_renorm real (r8), dimension(:,:,:,:), allocatable :: & TEMP_DATA ! temp array for reading monthly data type (datafile) :: & forcing_file ! io data file type for input forcing file type (io_field_desc) :: & io_tau, &! io field descriptor for input ws (both components) io_taux, &! io field descriptor for input ws in x direction io_tauy ! io field descriptor for input ws in y direction type (io_dim) :: & i_dim, j_dim, &! dimension descriptors for horiz dims month_dim ! dimension descriptor for monthly data !----------------------------------------------------------------------- ! ! read wind stress namelist input after setting default values. ! !----------------------------------------------------------------------- ws_data_type = 'analytic' ws_data_inc = 1.e20_r8 ws_interp_type = 'nearest' ws_interp_freq = 'never' ws_interp_inc = 1.e20_r8 ws_filename = 'unknown-ws' ws_file_fmt = 'bin' ws_data_renorm = c1 if (my_task == master_task) then open (nml_in, file=nml_filename, status='old',iostat=nml_error) if (nml_error /= 0) then nml_error = -1 else nml_error = 1 endif do while (nml_error > 0) read(nml_in, nml=forcing_ws_nml,iostat=nml_error) end do if (nml_error == 0) close(nml_in) endif call broadcast_scalar(nml_error, master_task) if (nml_error /= 0) then call exit_POP(sigAbort,'ERROR reading forcing_ws_nml') endif call broadcast_scalar(ws_data_type, master_task) call broadcast_scalar(ws_data_inc, master_task) call broadcast_scalar(ws_interp_type, master_task) call broadcast_scalar(ws_interp_freq, master_task) call broadcast_scalar(ws_interp_inc, master_task) call broadcast_scalar(ws_filename, master_task) call broadcast_scalar(ws_file_fmt, master_task) call broadcast_array (ws_data_renorm, master_task) ws_formulation = char_blank !----------------------------------------------------------------------- ! ! convert data_type to 'monthly-calendar' if input is 'monthly' ! !----------------------------------------------------------------------- if (ws_data_type == 'monthly') ws_data_type = 'monthly-calendar' !----------------------------------------------------------------------- ! ! convert interp_type to corresponding integer value. ! !----------------------------------------------------------------------- select case (ws_interp_type) case ('nearest') ws_interp_order = 1 case ('linear') ws_interp_order = 2 case ('4point') ws_interp_order = 4 case default call exit_POP(sigAbort, & 'init_ws: Unknown value for ws_interp_type') end select !----------------------------------------------------------------------- ! ! set values of the wind stress arrays (SMF or SMF_DATA) ! depending on the type of the wind stress data. ! !----------------------------------------------------------------------- select case (ws_data_type) !----------------------------------------------------------------------- ! ! no wind stress, therefore no interpolation in time is needed, ! nor are there any new values to be used. ! !----------------------------------------------------------------------- case ('none') SMF = c0 SMFT = c0 lsmft_avail = .true. ws_data_next = never ws_data_update = never ws_interp_freq = 'never' !----------------------------------------------------------------------- ! ! simple analytic wind stress that is constant in time, therefore no ! interpolation in time is needed, nor are there any new values ! to be used. ! !----------------------------------------------------------------------- case ('analytic') ws_data_next = never ws_data_update = never ws_interp_freq = 'never' !$OMP PARALLEL DO PRIVATE(iblock) do iblock = 1,nblocks_clinic ! simple zonal windstress field SMF (:,:,1,iblock) = -cos(3.0_r8*ULAT(:,:,iblock)) SMFT(:,:,1,iblock) = -cos(3.0_r8*TLAT(:,:,iblock)) ! Zero the zonal windstress at points within 1/100th degree ! of the true North Pole. ! NOTE this is needed because in tripole grids the ! true North Pole lies exactly on a U-point, and at ! that point ANGLE is not defined. ! In principle, this could also occur in a dipole grid ! if a U-point or T-point is near or very near the North Pole. ! This problem should also be dealt with in forcing_coupled.F ! and any other routines that use ANGLE and ANGLET at the ! North Pole point. !where(abs(abs(radian*ULAT(:,:,iblock)) - 90.0_r8) < 0.01_r8) & ! SMF(:,:,1,iblock) = c0 !where(abs(abs(radian*TLAT(:,:,iblock)) - 90.0_r8) < 0.01_r8) & ! SMFT(:,:,1,iblock) = c0 ! rotate vectors to grid SMF (:,:,2,iblock) = -sin(ANGLE (:,:,iblock))*SMF (:,:,1,iblock) SMFT(:,:,2,iblock) = -sin(ANGLET(:,:,iblock))*SMFT(:,:,1,iblock) SMF (:,:,1,iblock) = cos(ANGLE (:,:,iblock))*SMF (:,:,1,iblock) SMFT(:,:,1,iblock) = cos(ANGLET(:,:,iblock))*SMFT(:,:,1,iblock) if (ws_data_renorm(1) /= c1) then SMF (:,:,:,iblock) = SMF (:,:,:,iblock)*ws_data_renorm(1) SMFT(:,:,:,iblock) = SMFT(:,:,:,iblock)*ws_data_renorm(1) endif end do !$OMP END PARALLEL DO lsmft_avail = .true. !----------------------------------------------------------------------- ! ! annual mean climatological wind stress (read in from a file) that ! is constant in time, therefore no interpolation in time is needed, ! nor are there any new values to be used. ! !----------------------------------------------------------------------- case ('annual') ws_data_next = never ws_data_update = never ws_interp_freq = 'never' lsmft_avail = .false. allocate(TEMP_DATA(nx_block,ny_block,max_blocks_clinic,2)) allocate(ws_data_names(2), & ws_bndy_loc (2), & ws_bndy_type (2)) ws_data_names(1) = 'TAUX' ws_bndy_loc (1) = field_loc_NEcorner ws_bndy_type (1) = field_type_vector ws_data_names(2) = 'TAUY' ws_bndy_loc (2) = field_loc_NEcorner ws_bndy_type (2) = field_type_vector forcing_file = construct_file(ws_file_fmt, & full_name=trim(ws_filename), & record_length = rec_type_dbl, & recl_words=nx_global*ny_global) call data_set(forcing_file,'open_read') i_dim = construct_io_dim('i',nx_global) j_dim = construct_io_dim('j',ny_global) io_taux = construct_io_field(trim(ws_data_names(1)), & dim1=i_dim, dim2=j_dim, & field_loc = ws_bndy_loc(1), & field_type = ws_bndy_type(1), & d2d_array=TEMP_DATA(:,:,:,1)) io_tauy = construct_io_field(trim(ws_data_names(2)), & dim1=i_dim, dim2=j_dim, & field_loc = ws_bndy_loc(2), & field_type = ws_bndy_type(2), & d2d_array=TEMP_DATA(:,:,:,2)) call data_set(forcing_file,'define',io_taux) call data_set(forcing_file,'define',io_tauy) call data_set(forcing_file,'read' ,io_taux) call data_set(forcing_file,'read' ,io_tauy) call destroy_io_field(io_taux) call destroy_io_field(io_tauy) call data_set(forcing_file,'close') call destroy_file(forcing_file) !$OMP PARALLEL DO do iblock=1,nblocks_clinic SMF(:,:,1,iblock) = TEMP_DATA(:,:,iblock,1) SMF(:,:,2,iblock) = TEMP_DATA(:,:,iblock,2) if (ws_data_renorm(1) /= c1) SMF(:,:,:,iblock) = & SMF(:,:,:,iblock)*ws_data_renorm(1) end do !$OMP END PARALLEL DO deallocate(TEMP_DATA) if (my_task == master_task) then write(stdout,blank_fmt) write(stdout,'(a22,a)') ' WS Annual file read: ', & trim(ws_filename) endif !----------------------------------------------------------------------- ! ! monthly mean climatological wind stress. all 12 months are read ! in from a file. interpolation order (ws_interp_order) may be ! specified with namelist input. ! !----------------------------------------------------------------------- case ('monthly-equal','monthly-calendar') lsmft_avail = .false. allocate( SMF_DATA(nx_block,ny_block,max_blocks_clinic,2,0:12), & TEMP_DATA(nx_block,ny_block,12,max_blocks_clinic)) allocate(ws_data_names(2), & ws_bndy_loc (2), & ws_bndy_type (2)) SMF_DATA = c0 ws_data_names(1) = 'TAUX' ws_bndy_loc (1) = field_loc_NEcorner ws_bndy_type (1) = field_type_vector ws_data_names(2) = 'TAUY' ws_bndy_loc (2) = field_loc_NEcorner ws_bndy_type (2) = field_type_vector call find_forcing_times(ws_data_time, ws_data_inc, & ws_interp_type, ws_data_next, & ws_data_time_min_loc, ws_data_update, & ws_data_type) forcing_file = construct_file(ws_file_fmt, & full_name=trim(ws_filename), & record_length = rec_type_dbl, & recl_words=nx_global*ny_global) call data_set(forcing_file,'open_read') i_dim = construct_io_dim('i',nx_global) j_dim = construct_io_dim('j',ny_global) month_dim = construct_io_dim('month',12) io_taux = construct_io_field(trim(ws_data_names(1)), & dim1=i_dim, dim2=j_dim, dim3=month_dim, & field_loc = ws_bndy_loc(1), & field_type = ws_bndy_type(1), & d3d_array=TEMP_DATA) io_tauy = construct_io_field(trim(ws_data_names(2)), & dim1=i_dim, dim2=j_dim, dim3=month_dim, & field_loc = ws_bndy_loc(2), & field_type = ws_bndy_type(2), & d3d_array=TEMP_DATA) call data_set(forcing_file,'define',io_taux) call data_set(forcing_file,'define',io_tauy) call data_set(forcing_file,'read',io_taux) !$OMP PARALLEL DO PRIVATE(iblock, n) do iblock=1,nblocks_clinic do n=1,12 SMF_DATA(:,:,iblock,1,n) = TEMP_DATA(:,:,n,iblock) end do end do !$OMP END PARALLEL DO call destroy_io_field(io_taux) call data_set(forcing_file,'read',io_tauy) !$OMP PARALLEL DO PRIVATE(iblock, n) do iblock=1,nblocks_clinic do n=1,12 SMF_DATA(:,:,iblock,2,n) = TEMP_DATA(:,:,n,iblock) end do end do !$OMP END PARALLEL DO call destroy_io_field(io_tauy) call data_set(forcing_file,'close') call destroy_file(forcing_file) deallocate(TEMP_DATA) if (ws_data_renorm(1) /= c1) SMF_DATA = SMF_DATA*ws_data_renorm(1) if (my_task == master_task) then write(stdout,blank_fmt) write(stdout,'(a23,a)') ' WS Monthly file read: ', & trim(ws_filename) endif !----------------------------------------------------------------------- ! ! wind stress specified every n-hours, where the n-hour increment ! should be specified with namelist input (ws_data_inc). ! only as many times as are necessary based on the order ! of the temporal interpolation scheme (ws_interp_order) reside ! in memory at any given time. ! !----------------------------------------------------------------------- case ('n-hour') lsmft_avail = .false. allocate( SMF_DATA(nx_block,ny_block,max_blocks_clinic,2, & 0:ws_interp_order)) allocate(ws_data_names(2), & ws_bndy_loc (2), & ws_bndy_type (2)) SMF_DATA = c0 ws_data_names(1) = 'TAUX' ws_bndy_loc (1) = field_loc_NEcorner ws_bndy_type (1) = field_type_vector ws_data_names(2) = 'TAUY' ws_bndy_loc (2) = field_loc_NEcorner ws_bndy_type (2) = field_type_vector call find_forcing_times(ws_data_time, ws_data_inc, & ws_interp_type, ws_data_next, & ws_data_time_min_loc, ws_data_update, & ws_data_type) do n = 1, ws_interp_order call get_forcing_filename(forcing_filename, ws_filename, & ws_data_time(n), ws_data_inc) forcing_file = construct_file(ws_file_fmt, & full_name=trim(forcing_filename), & record_length = rec_type_dbl, & recl_words=nx_global*ny_global) call data_set(forcing_file,'open_read') i_dim = construct_io_dim('i',nx_global) j_dim = construct_io_dim('j',ny_global) io_taux = construct_io_field(trim(ws_data_names(1)), & dim1=i_dim, dim2=j_dim, & field_loc = ws_bndy_loc(1), & field_type = ws_bndy_type(1), & d2d_array=SMF_DATA(:,:,:,1,n)) io_tauy = construct_io_field(trim(ws_data_names(2)), & dim1=i_dim, dim2=j_dim, & field_loc = ws_bndy_loc(2), & field_type = ws_bndy_type(2), & d2d_array=SMF_DATA(:,:,:,2,n)) call data_set(forcing_file,'define',io_taux) call data_set(forcing_file,'define',io_tauy) call data_set(forcing_file,'read' ,io_taux) call data_set(forcing_file,'read' ,io_tauy) call destroy_io_field(io_taux) call destroy_io_field(io_tauy) call data_set(forcing_file,'close') call destroy_file(forcing_file) if (my_task == master_task) then write(stdout,blank_fmt) write(stdout,'(a22,a)') ' WS n-hour file read: ', & trim(forcing_filename) endif enddo if (ws_data_renorm(1) /= c1) SMF_DATA = SMF_DATA*ws_data_renorm(1) case default call exit_POP(sigAbort,'init_ws: Unknown value for ws_data_type') end select !----------------------------------------------------------------------- ! ! now check interpolation period (ws_interp_freq) to set the ! time for the next temporal interpolation (ws_interp_next). ! ! if no interpolation is to be done, set next interpolation time ! to a large number so the wind stress update test in routine ! set_surface_forcing will always be false. ! ! if interpolation is to be done every n-hours, find the first ! interpolation time greater than the current time. ! ! if interpolation is to be done every timestep, set next interpolation ! time to a large negative number so the wind stress update test in ! routine set_surface_forcing will always be true. ! !----------------------------------------------------------------------- select case (ws_interp_freq) case ('never') ws_interp_next = never ws_interp_last = never ws_interp_inc = c0 case ('n-hour') call find_interp_time(ws_interp_inc, ws_interp_next) case ('every-timestep') ws_interp_next = always ws_interp_inc = c0 case default call exit_POP(sigAbort, & 'init_ws: Unknown value for ws_interp_freq') end select if (nsteps_total == 0) ws_interp_last = thour00 !----------------------------------------------------------------------- ! ! echo forcing options to stdout. ! !----------------------------------------------------------------------- ws_data_label = 'Surface Wind Stress' call echo_forcing_options(ws_data_type, & ws_formulation, ws_data_inc, & ws_interp_freq, ws_interp_type, & ws_interp_inc, ws_data_label) !----------------------------------------------------------------------- !EOC end subroutine init_ws !*********************************************************************** !BOP ! !IROUTINE: set_ws ! !INTERFACE: subroutine set_ws(SMF,SMFT) ! !DESCRIPTION: ! Update the current value of the wind stress arrays (SMF) by ! interpolating to the current time. It may also be necessary to ! use new data values than were used in the previous interpolation. ! ! !REVISION HISTORY: ! same as module ! !INPUT/OUTPUT PARAMETERS: real (r8), dimension(nx_block,ny_block,2,max_blocks_clinic), & intent(inout) :: & SMF ! surface momentum fluxes at current timestep real (r8), dimension(nx_block,ny_block,2,max_blocks_clinic), & intent(inout), optional :: & SMFT ! surface momentum fluxes at T points if available !EOP !BOC !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- integer (int_kind) :: & iblock, & first, second, third, fourth, & ws_data_time_max_loc !----------------------------------------------------------------------- ! ! check if new data is necessary for interpolation. if yes, then ! shuffle indices in SMF_DATA, and ws_data_time arrays ! and read in new data if necessary ('n-hour' case). also ! increment values of ws_data_time_min_loc, ws_data_next and ! ws_data_update.then do the temporal interpolation. ! !----------------------------------------------------------------------- if (thour00 >= ws_interp_next .or. nsteps_run == 0) then select case(ws_data_type) case ('monthly-equal','monthly-calendar') ws_data_label = 'WS Monthly' if (thour00 >= ws_data_update) & call update_forcing_data(ws_data_time, ws_data_time_min_loc,& ws_interp_type, ws_data_next, & ws_data_update, ws_data_type, & ws_data_inc, SMF_DATA(:,:,:,:,1:12), & ws_data_renorm, & ws_data_label, ws_data_names, & ws_bndy_loc, ws_bndy_type, & ws_filename, ws_file_fmt) call interpolate_forcing(SMF_DATA(:,:,:,:,0), & SMF_DATA(:,:,:,:,1:12), & ws_data_time, ws_interp_type, & ws_data_time_min_loc, ws_interp_freq, & ws_interp_inc, ws_interp_next, & ws_interp_last, nsteps_run) !$OMP PARALLEL DO do iblock=1,nblocks_clinic SMF(:,:,1,iblock) = SMF_DATA(:,:,iblock,1,0) SMF(:,:,2,iblock) = SMF_DATA(:,:,iblock,2,0) end do !$OMP END PARALLEL DO case ('n-hour') ws_data_label = 'WS n-hour' if (thour00 >= ws_data_update) & call update_forcing_data(ws_data_time, ws_data_time_min_loc,& ws_interp_type, ws_data_next, & ws_data_update, ws_data_type, & ws_data_inc, & SMF_DATA(:,:,:,:,1:ws_interp_order), & ws_data_renorm, & ws_data_label, ws_data_names, & ws_bndy_loc, ws_bndy_type, & ws_filename, ws_file_fmt) call interpolate_forcing(SMF_DATA(:,:,:,:,0), & SMF_DATA(:,:,:,:,1:ws_interp_order), & ws_data_time, ws_interp_type, & ws_data_time_min_loc, ws_interp_freq, & ws_interp_inc, ws_interp_next, & ws_interp_last, nsteps_run) !$OMP PARALLEL DO do iblock=1,nblocks_clinic SMF(:,:,1,iblock) = SMF_DATA(:,:,iblock,1,0) SMF(:,:,2,iblock) = SMF_DATA(:,:,iblock,2,0) end do !$OMP END PARALLEL DO end select if (nsteps_run /= 0) ws_interp_next = ws_interp_next + & ws_interp_inc endif ! thour00 > interp_next !----------------------------------------------------------------------- !EOC end subroutine set_ws !*********************************************************************** end module forcing_ws !|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||