!||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| module forcing_coupled !BOP !MODULE: forcing_coupled ! !DESCRIPTION: ! This module contains all the routines necessary for coupling POP to ! atmosphere and sea ice models using the NCAR CCSM flux coupler. To ! enable the routines in this module, the coupled ifdef option must ! be specified during the make process. ! ! !REVISION HISTORY: ! SVN:$Id: forcing_coupled.F90 77949 2016-04-13 16:26:53Z jet $ ! ! !USES: use POP_KindsMod use POP_ErrorMod use POP_CommMod use POP_FieldMod use POP_GridHorzMod use POP_HaloMod use kinds_mod use blocks, only: nx_block, ny_block, block, get_block use domain_size use domain use io_types, only: stdout, nml_in use communicate use global_reductions use constants use io use time_management use grid use prognostic use exit_mod use ice, only: tfreez, tmelt, liceform,QFLUX, QICE, AQICE, tlast_ice use forcing_shf use forcing_sfwf use forcing_ws, only: ws_data_type use forcing_fields use timers !*** ccsm use ms_balance use tavg use registry use named_field_mod, only: named_field_register, named_field_get_index, & named_field_set, named_field_get use forcing_fields use estuary_mod use mcog, only: tavg_mcog use mcog, only: FRAC_BIN, QSW_RAW_BIN implicit none save !EOP !BOC !----------------------------------------------------------------------- ! ! module variables ! !----------------------------------------------------------------------- integer (int_kind) :: & coupled_freq_iopt, &! coupler frequency option coupled_freq, &! frequency of coupling ncouple_per_day ! num of coupler comms per day #if CCSMCOUPLED !----------------------------------------------------------------------- ! ! ids for tavg diagnostics computed from forcing_coupled ! !----------------------------------------------------------------------- integer (int_kind) :: & tavg_EVAP_F, &! tavg id for evaporation flux tavg_PREC_F, &! tavg id for precipitation flux (rain + snow) tavg_SNOW_F, &! tavg id for snow flux tavg_MELT_F, &! tavg id for melt flux tavg_ROFF_F, &! tavg id for river runoff flux tavg_IOFF_F, &! tavg id for ice runoff flux due to land-model snow capping tavg_SALT_F, &! tavg id for salt flux tavg_SENH_F, &! tavg id for sensible heat flux tavg_LWUP_F, &! tavg id for longwave heat flux up tavg_LWDN_F, &! tavg id for longwave heat flux dn tavg_MELTH_F, &! tavg id for melt heat flux tavg_IFRAC ! tavg id for ice fraction #endif !----------------------------------------------------------------------- ! ! Options for distributing net shortwave heat flux over a coupling ! interval. All options preserve time-integrated flux. ! !----------------------------------------------------------------------- integer (int_kind), parameter :: & qsw_distrb_iopt_const = 1, &! qsw constant over a coupling interval qsw_distrb_iopt_12hr = 2, &! qsw smoothly spread over 12 hour window ! only works for daily coupling qsw_distrb_iopt_cosz = 3 ! qsw proportional to cos of solar zenith angle integer (int_kind) :: qsw_distrb_iopt real (r8), dimension(:), allocatable :: & qsw_12hr_factor !----------------------------------------------------------------------- ! variables for qsw cosz option !----------------------------------------------------------------------- integer (int_kind) :: timer_compute_cosz real (r8) :: & tday00_interval_beg, & ! model time at beginning of coupling interval orb_eccen, & ! Earth eccentricity orb_obliqr, & ! Earth Obliquity orb_lambm0, & ! longitude of perihelion at v-equinox orb_mvelpp ! Earths Moving vernal equinox of orbit +pi real (r8), dimension(:,:,:), allocatable :: & QSW_COSZ_WGHT, & ! weights QSW_COSZ_WGHT_NORM ! normalization for QSW_COSZ_WGHT integer (int_kind), private :: & cpl_ts ! flag id for coupled_ts flag !EOC !*********************************************************************** contains !*********************************************************************** !BOP ! !IROUTINE: pop_init_coupled ! !INTERFACE: subroutine pop_init_coupled ! !DESCRIPTION: ! This routine sets up everything necessary for coupling with CCSM4. ! ! !REVISION HISTORY: ! same as module !EOP !BOC !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- character (char_len) :: & coupled_freq_opt, qsw_distrb_opt namelist /coupled_nml/ coupled_freq_opt, coupled_freq, qsw_distrb_opt integer (int_kind) :: & k, iblock, nsend, & nml_error ! namelist i/o error flag type (block) :: & this_block ! block information for current block !----------------------------------------------------------------------- ! ! variables associated with qsw 12hr ! !----------------------------------------------------------------------- real (r8) :: & time_for_forcing, &! time of day for surface forcing frac_day_forcing, &! fraction of day based on time_for_forcing cycle_function, &! intermediate result weight_forcing, &! forcing weights sum_forcing ! sum of forcing weights integer (int_kind) :: & count_forcing ! time step counter (== nsteps_this_interval+1) integer (int_kind) :: & i,j,n !----------------------------------------------------------------------- ! ! read coupled_nml namelist to start coupling and determine ! coupling frequency ! !----------------------------------------------------------------------- coupled_freq_opt = 'never' coupled_freq_iopt = freq_opt_never coupled_freq = 100000 qsw_distrb_opt = 'const' 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=coupled_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 coupled_nml') endif if (my_task == master_task) then write(stdout,blank_fmt) write(stdout,ndelim_fmt) write(stdout,blank_fmt) write(stdout,*) ' Coupling:' write(stdout,blank_fmt) write(stdout,*) ' coupled_nml namelist settings:' write(stdout,blank_fmt) write(stdout, coupled_nml) write(stdout,blank_fmt) endif if (my_task == master_task) then select case (coupled_freq_opt) case ('nyear') coupled_freq_iopt = -1000 case ('nmonth') coupled_freq_iopt = -1000 case ('nday') if (coupled_freq == 1) then coupled_freq_iopt = freq_opt_nday ncouple_per_day = 1 else coupled_freq_iopt = -1000 endif case ('nhour') if (coupled_freq <= 24) then coupled_freq_iopt = freq_opt_nhour ncouple_per_day = 24/coupled_freq else coupled_freq_iopt = -1000 endif case ('nsecond') if (coupled_freq <= seconds_in_day) then coupled_freq_iopt = freq_opt_nsecond ncouple_per_day = seconds_in_day/coupled_freq else coupled_freq_iopt = -1000 endif case ('nstep') if (coupled_freq <= nsteps_per_day) then coupled_freq_iopt = freq_opt_nstep ncouple_per_day = nsteps_per_day/coupled_freq else coupled_freq_iopt = -1000 endif case ('never') coupled_freq_iopt = -9999 case default coupled_freq_iopt = -2000 end select select case (qsw_distrb_opt) case ('const') qsw_distrb_iopt = qsw_distrb_iopt_const case ('12hr') qsw_distrb_iopt = qsw_distrb_iopt_12hr case ('cosz') qsw_distrb_iopt = qsw_distrb_iopt_cosz call register_string('qsw_distrb_iopt_cosz') case default qsw_distrb_iopt = -1000 end select endif call broadcast_scalar(coupled_freq_iopt, master_task) call broadcast_scalar(coupled_freq , master_task) call broadcast_scalar(qsw_distrb_iopt , master_task) call broadcast_scalar(ncouple_per_day , master_task) if (coupled_freq_iopt == -1000) then call exit_POP(sigAbort, & 'ERROR: Coupling frequency must be at least once per day') else if (coupled_freq_iopt == -2000) then call exit_POP(sigAbort, & 'ERROR: Unknown option for coupling frequency') endif if (registry_match('lcoupled') .eqv. (coupled_freq_iopt == -9999) ) then call exit_POP(sigAbort, & 'ERROR: inconsistency between lcoupled and coupled_freq_iopt settings') endif if (qsw_distrb_iopt == -1000) then call exit_POP(sigAbort, & 'ERROR: Unknown option for qsw_distrb_opt') endif !----------------------------------------------------------------------- ! ! check consistency of the qsw_distrb_iopt option with various ! time manager options ! !----------------------------------------------------------------------- if ( (qsw_distrb_iopt == qsw_distrb_iopt_12hr) .or. & (qsw_distrb_iopt == qsw_distrb_iopt_cosz) ) then if ( tmix_iopt == tmix_avgfit .or. tmix_iopt == tmix_robert) then ! ok; these options are supported else call exit_POP(sigAbort, & 'ERROR: time_mix_opt must be set to avgfit for qsw_distrb_opt '/& &/ 'of 12hr or cosz') endif if ( dttxcel(1) /= c1 .or. dtuxcel /= c1 ) & call exit_POP(sigAbort, & 'ERROR: using the specified accelerated integration '/& &/ 'technique may not be appropriate for qsw_distrb_opt '/& &/ 'of 12hr or cosz') endif !----------------------------------------------------------------------- ! ! allocate and compute the short wave heat flux multiplier for qsw 12hr ! !----------------------------------------------------------------------- allocate ( qsw_12hr_factor(nsteps_per_interval)) qsw_12hr_factor = c1 if ( qsw_distrb_iopt == qsw_distrb_iopt_12hr ) then ! mimic a day time_for_forcing = c0 count_forcing = 1 sum_forcing = c0 do n=1,nsteps_per_interval frac_day_forcing = time_for_forcing / seconds_in_day cycle_function = cos( pi * ( c2 * frac_day_forcing - c1 ) ) qsw_12hr_factor(n) = c2 * ( cycle_function & + abs(cycle_function) ) & * cycle_function weight_forcing = c1 if ( count_forcing == 2 .or. & mod(count_forcing,time_mix_freq) == 0 ) & weight_forcing = p5 time_for_forcing = time_for_forcing + weight_forcing * dt(1) sum_forcing = sum_forcing & + weight_forcing * dt(1) * qsw_12hr_factor(n) count_forcing = count_forcing + 1 enddo qsw_12hr_factor = qsw_12hr_factor * seconds_in_day & / sum_forcing ! check the final integral count_forcing = 1 sum_forcing = c0 do n=1,nsteps_per_interval weight_forcing = c1 if ( count_forcing == 2 .or. & mod(count_forcing,time_mix_freq) == 0 ) & weight_forcing = p5 sum_forcing = sum_forcing & + weight_forcing * dt(1) * qsw_12hr_factor(n) count_forcing = count_forcing + 1 enddo if ( sum_forcing < (seconds_in_day - 1.0e-5_r8) .or. & sum_forcing > (seconds_in_day + 1.0e-5_r8) ) & call exit_POP (sigAbort, & 'ERROR: qsw 12hr temporal integral is incorrect') endif !----------------------------------------------------------------------- ! ! allocate space for qsw cosz fields ! !----------------------------------------------------------------------- if ( qsw_distrb_iopt == qsw_distrb_iopt_cosz ) then allocate( & QSW_COSZ_WGHT(nx_block,ny_block,nblocks_clinic), & QSW_COSZ_WGHT_NORM(nx_block,ny_block,nblocks_clinic)) endif #if CCSMCOUPLED !----------------------------------------------------------------------- ! ! define tavg fields computed from forcing_coupled routines ! !----------------------------------------------------------------------- call define_tavg_field(tavg_EVAP_F,'EVAP_F',2, & long_name='Evaporation Flux from Coupler', & units='kg/m^2/s', grid_loc='2110', & coordinates='TLONG TLAT time') call define_tavg_field(tavg_PREC_F,'PREC_F',2, & long_name='Precipitation Flux from Cpl (rain+snow)', & units='kg/m^2/s', grid_loc='2110', & coordinates='TLONG TLAT time') call define_tavg_field(tavg_SNOW_F,'SNOW_F',2, & long_name='Snow Flux from Coupler', & units='kg/m^2/s', grid_loc='2110', & coordinates='TLONG TLAT time') call define_tavg_field(tavg_MELT_F,'MELT_F',2, & long_name='Melt Flux from Coupler', & units='kg/m^2/s', grid_loc='2110', & coordinates='TLONG TLAT time') call define_tavg_field(tavg_ROFF_F,'ROFF_F',2, & long_name='Runoff Flux from Coupler', & units='kg/m^2/s', grid_loc='2110', & coordinates='TLONG TLAT time') call define_tavg_field(tavg_IOFF_F,'IOFF_F',2, & long_name='Ice Runoff Flux from Coupler due to Land-Model Snow Capping', & units='kg/m^2/s', grid_loc='2110', & coordinates='TLONG TLAT time') call define_tavg_field(tavg_SALT_F,'SALT_F',2, & long_name='Salt Flux from Coupler (kg of salt/m^2/s)',& units='kg/m^2/s', grid_loc='2110', & coordinates='TLONG TLAT time') call define_tavg_field(tavg_SENH_F,'SENH_F',2, & long_name='Sensible Heat Flux from Coupler', & units='watt/m^2', grid_loc='2110', & coordinates='TLONG TLAT time') call define_tavg_field(tavg_LWUP_F,'LWUP_F',2, & long_name='Longwave Heat Flux (up) from Coupler', & units='watt/m^2', grid_loc='2110', & coordinates='TLONG TLAT time') call define_tavg_field(tavg_LWDN_F,'LWDN_F',2, & long_name='Longwave Heat Flux (dn) from Coupler', & units='watt/m^2', grid_loc='2110', & coordinates='TLONG TLAT time') call define_tavg_field(tavg_MELTH_F,'MELTH_F',2, & long_name='Melt Heat Flux from Coupler', & units='watt/m^2', grid_loc='2110', & coordinates='TLONG TLAT time') call define_tavg_field(tavg_IFRAC,'IFRAC',2, & long_name='Ice Fraction from Coupler', & units='fraction', grid_loc='2110', & coordinates='TLONG TLAT time') !----------------------------------------------------------------------- ! ! Initialize flags and shortwave absorption profile ! Note that the cpl_write_xxx flags have _no_ default value; ! therefore, they must be explicitly set .true. and .false. ! at the appropriate times ! !----------------------------------------------------------------------- call init_time_flag('coupled_ts', cpl_ts, & owner='pop_init_coupled', & freq_opt = coupled_freq_iopt, & freq = coupled_freq) !----------------------------------------------------------------------- ! ! If this is a restart, then read_restart knows the last timestep was ! a coupled timestep and has registered the string 'coupled_ts_last_true' ! (read_restart was called prior to the initialization of coupled_ts) ! !----------------------------------------------------------------------- if (registry_match('coupled_ts_last_true') ) & call override_time_flag (cpl_ts, old_value=.true.) lsmft_avail = .true. !----------------------------------------------------------------------- ! ! initialize timer for computing cosz ! !----------------------------------------------------------------------- if ( qsw_distrb_iopt == qsw_distrb_iopt_cosz ) then call get_timer (timer_compute_cosz, 'COMPUTE_COSZ', nblocks_clinic, & distrb_clinic%nprocs) endif !----------------------------------------------------------------------- ! ! register this subroutine ! !----------------------------------------------------------------------- call register_string('pop_init_coupled') #endif !----------------------------------------------------------------------- !EOC call flushm (stdout) end subroutine pop_init_coupled !*********************************************************************** !BOP ! !IROUTINE: pop_init_partially_coupled ! !INTERFACE: subroutine pop_init_partially_coupled ! !DESCRIPTION: ! This routine initializes and allocates arrays for the partially-coupled ! option ! ! !REVISION HISTORY: ! same as module !EOP !BOC #if CCSMCOUPLED !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- logical (log_kind) :: & lcoupled character (char_len) :: & message integer (int_kind) :: & number_of_fatal_errors lcoupled = registry_match('lcoupled') if ( lcoupled .and. shf_formulation /= 'partially-coupled' ) then shf_num_comps = 1 shf_comp_qsw = 1 allocate(SHF_COMP(nx_block,ny_block,max_blocks_clinic,shf_num_comps)) SHF_COMP = c0 endif !----------------------------------------------------------------------- ! ! initialize and allocate some partially coupled variables ! !----------------------------------------------------------------------- if ( lcoupled & .and. sfwf_formulation /= 'partially-coupled' & .and. sfc_layer_type == sfc_layer_varthick .and. & .not. lfw_as_salt_flx .and. liceform ) then sfwf_num_comps = 1 sfwf_comp_cpl = 1 tfw_num_comps = 1 tfw_comp_cpl = 1 allocate(SFWF_COMP(nx_block,ny_block, max_blocks_clinic,sfwf_num_comps)) allocate( TFW_COMP(nx_block,ny_block,nt,max_blocks_clinic, tfw_num_comps)) SFWF_COMP = c0 TFW_COMP = c0 endif !----------------------------------------------------------------------- ! ! check compatibility of partially-coupled option with other options ! !----------------------------------------------------------------------- number_of_fatal_errors = 0 if (.not. lcoupled .and. (shf_formulation == 'partially-coupled' .or. & sfwf_formulation == 'partially-coupled' ) ) then message = & 'ERROR: partially-coupled option is allowed only when coupled' write(stdout,*) message number_of_fatal_errors = number_of_fatal_errors + 1 endif if (lcoupled .and. (shf_formulation == 'partially-coupled' .and. & sfwf_formulation /= 'partially-coupled') .or. & (shf_formulation /= 'partially-coupled' .and. & sfwf_formulation == 'partially-coupled') ) then message = & 'partially-coupled must be used for both shf and sfwf' write(stdout,*) message number_of_fatal_errors = number_of_fatal_errors + 1 endif if (lcoupled .and. shf_formulation /= 'partially-coupled' .and. & shf_data_type /= 'none') then message = & 'shf_data_type must be set to none or '/& &/ 'shf_formulation must be partially_coupled when lcoupled is true' write(stdout,*) message number_of_fatal_errors = number_of_fatal_errors + 1 endif if (lcoupled .and. sfwf_formulation /= 'partially-coupled' .and. & sfwf_data_type /= 'none') then message = & 'sfwf_data_type must be set to none or '/& &/ 'sfwf_formulation must be partially_coupled when lcoupled is true' write(stdout,*) message number_of_fatal_errors = number_of_fatal_errors + 1 endif !----------------------------------------------------------------------- ! ! check coupled compatibility with other forcing options ! !----------------------------------------------------------------------- if (lcoupled .and. ws_data_type /= 'none') then message = & 'ws_data_type must be set to none in coupled mode' write(stdout,*) message number_of_fatal_errors = number_of_fatal_errors + 1 endif if (number_of_fatal_errors /= 0) & call exit_POP(sigAbort,'subroutine pop_init_partially_coupled') #endif !----------------------------------------------------------------------- !EOC call flushm (stdout) end subroutine pop_init_partially_coupled !*********************************************************************** !BOP ! !IROUTINE: pop_set_coupled_forcing ! !INTERFACE: subroutine pop_set_coupled_forcing ! !DESCRIPTION: ! This routine is called immediately following the receipt of fluxes ! from the coupler. It combines fluxes received from the coupler into ! the STF array and converts from W/m**2 into model units. It also ! balances salt/freshwater in marginal seas and sets SHF_QSW_RAW ! and SHF_COMP. Compute QSW_COSZ_WGHT_NORM if needed. ! ! !REVISION HISTORY: ! same as module !EOP !BOC !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- #if CCSMCOUPLED integer (int_kind) :: n, nn, iblock real (r8) :: cosz_day ! time where cosz is computed real (r8), dimension(nx_block,ny_block,max_blocks_clinic) :: & WORK1, WORK2 ! local work space !----------------------------------------------------------------------- ! ! combine heat flux components into STF array and convert from W/m**2 ! (note: latent heat flux = evaporation*latent_heat_vapor_mks) ! (note: snow melt heat flux = - snow_f*latent_heat_fusion_mks) ! !----------------------------------------------------------------------- !*** need to zero out any padded cells WORK1 = c0 WORK2 = c0 !$OMP PARALLEL DO PRIVATE(iblock) do iblock = 1, nblocks_clinic STF(:,:,1,iblock) = (EVAP_F(:,:,iblock)*latent_heat_vapor_mks & + SENH_F(:,:,iblock) + LWUP_F(:,:,iblock) & + LWDN_F(:,:,iblock) + MELTH_F(:,:,iblock) & -(SNOW_F(:,:,iblock)+IOFF_F(:,:,iblock)) * latent_heat_fusion_mks)* & RCALCT(:,:,iblock)*hflux_factor enddo !$OMP END PARALLEL DO !----------------------------------------------------------------------- ! ! combine freshwater flux components ! ! for variable thickness surface layer, compute fresh water and ! salt fluxes ! !----------------------------------------------------------------------- if (sfc_layer_type == sfc_layer_varthick .and. & .not. lfw_as_salt_flx) then !*** compute fresh water flux (cm/s) !$OMP PARALLEL DO PRIVATE(iblock,n) do iblock = 1, nblocks_clinic FW(:,:,iblock) = RCALCT(:,:,iblock) * & ( PREC_F(:,:,iblock)+EVAP_F(:,:,iblock) & +ROFF_F(:,:,iblock)+IOFF_F(:,:,iblock))*fwmass_to_fwflux WORK1(:,:,iblock) = RCALCT(:,:,iblock) * & MELT_F(:,:,iblock) * fwmass_to_fwflux !*** compute tracer concentration in fresh water !*** in principle, temperature of each water flux !*** could be different. e.g. !TFW(:,:,1,iblock) = RCALCT(:,:,iblock)*fwmass_to_fwflux & ! (PREC_F(:,:,iblock)*TEMP_PREC(:,:,iblock) + & ! EVAP_F(:,:,iblock)*TEMP_EVAP(:,:,iblock) + & ! MELT_F(:,:,iblock)*TEMP_MELT(:,:,iblock) + & ! ROFF_F(:,:,iblock)*TEMP_ROFF(:,:,iblock)) !*** currently assume water comes in at sea surface temp call tmelt(WORK2(:,:,iblock),TRACER(:,:,1,2,curtime,iblock)) TFW(:,:,1,iblock) = FW(:,:,iblock)*TRACER(:,:,1,1,curtime,iblock) & + WORK1(:,:,iblock) * WORK2(:,:,iblock) FW(:,:,iblock) = FW(:,:,iblock) + WORK1(:,:,iblock) !*** compute salt flux !*** again, salinity could be different for each !*** component of water flux !TFW(:,:,2,iblock) = RCALCT(:,:,iblock)*fwmass_to_fwflux & ! (PREC_F(:,:,iblock)*SALT_PREC(:,:,iblock) + & ! EVAP_F(:,:,iblock)*SALT_EVAP(:,:,iblock) + & ! MELT_F(:,:,iblock)*SALT_MELT(:,:,iblock) + & ! ROFF_F(:,:,iblock)*SALT_ROFF(:,:,iblock)) !*** currently assume prec, evap and roff are fresh !*** and all salt come from ice melt where (MELT_F(:,:,iblock) /= c0) WORK1(:,:,iblock) = & SALT_F(:,:,iblock)/MELT_F(:,:,iblock) ! salinity (msu) of melt water elsewhere WORK1(:,:,iblock) = c0 end where TFW(:,:,2,iblock) = RCALCT(:,:,iblock)*MELT_F(:,:,iblock)* & fwmass_to_fwflux*WORK1(:,:,iblock) ! + PREC_F(:,:,iblock)*c0 + EVAP_F(:,:,iblock)*c0 + ROFF_F(:,:,iblock)*c0 + IOFF_F(:,:,iblock)*c0 do n=3,nt TFW(:,:,n,iblock) = c0 ! no additional tracers in fresh water end do enddo !$OMP END PARALLEL DO else ! convert fresh water to virtual salinity flux !----------------------------------------------------------------------- ! ! if not a variable thickness surface layer or if fw_as_salt_flx ! flag is on, convert fresh and salt inputs to a virtual salinity flux ! !----------------------------------------------------------------------- if ( lestuary_on) then ! Treat river runoff as the interior source if (lvsf_river) call set_estuary_vsf_forcing ! Include estuary exchange flow as vertical salt flux if (lebm_on) call set_estuary_exch_circ endif if ( lestuary_on ) then !$OMP PARALLEL DO PRIVATE(iblock) do iblock = 1, nblocks_clinic ! Add up the freshwater flux except river runoff STF(:,:,2,iblock) = RCALCT(:,:,iblock)*( & (PREC_F(:,:,iblock)+EVAP_F(:,:,iblock)+ & MELT_F(:,:,iblock)& +IOFF_F(:,:,iblock))*salinity_factor & + SALT_F(:,:,iblock)*sflux_factor) enddo !$OMP END PARALLEL DO if (lvsf_river) THEN !$OMP PARALLEL DO PRIVATE(iblock) ! Add global correction for salt conservation with local river reference salinity do iblock = 1, nblocks_clinic STF(:,:,2,iblock) = & & STF(:,:,2,iblock) + RCALCT(:,:,iblock)*vsf_river_correction enddo !$OMP END PARALLEL DO endif else !$OMP PARALLEL DO PRIVATE(iblock) ! Estuary exchange flow is not on do iblock = 1, nblocks_clinic STF(:,:,2,iblock) = RCALCT(:,:,iblock)*( & (PREC_F(:,:,iblock)+EVAP_F(:,:,iblock)+ & MELT_F(:,:,iblock)+ROFF_F(:,:,iblock)+IOFF_F(:,:,iblock))*salinity_factor & + SALT_F(:,:,iblock)*sflux_factor) enddo !$OMP END PARALLEL DO endif !----------------------------------------------------------------------- ! ! balance salt/freshwater in marginal seas ! !----------------------------------------------------------------------- if (lms_balance .and. sfwf_formulation /= 'partially-coupled' ) then call ms_balancing (STF(:,:,2,:),EVAP_F, PREC_F, MELT_F,ROFF_F,IOFF_F, & SALT_F, QFLUX, 'salt') endif endif !$OMP PARALLEL DO PRIVATE(iblock,n) do iblock = 1, nblocks_clinic SHF_QSW_RAW(:,:,iblock) = SHF_QSW(:,:,iblock) if ( shf_formulation == 'partially-coupled' ) then SHF_COMP(:,:,iblock,shf_comp_cpl) = STF(:,:,1,iblock) if ( .not. lms_balance ) then SHF_COMP(:,:,iblock,shf_comp_cpl) = & SHF_COMP(:,:,iblock,shf_comp_cpl) * MASK_SR(:,:,iblock) SHF_QSW(:,:,iblock) = SHF_QSW(:,:,iblock) * MASK_SR(:,:,iblock) endif endif SHF_COMP(:,:,iblock,shf_comp_qsw) = SHF_QSW(:,:,iblock) if ( sfwf_formulation == 'partially-coupled' ) then if (sfc_layer_type == sfc_layer_varthick .and. & .not. lfw_as_salt_flx) then SFWF_COMP(:,:,iblock,sfwf_comp_cpl) = & FW(:,:,iblock) * MASK_SR(:,:,iblock) do n=1,nt TFW_COMP(:,:,n,iblock,tfw_comp_cpl) = & TFW(:,:,n,iblock) * MASK_SR(:,:,iblock) enddo else SFWF_COMP(:,:,iblock,sfwf_comp_cpl) = & STF(:,:,2,iblock) * MASK_SR(:,:,iblock) endif else if ( sfc_layer_type == sfc_layer_varthick .and. & .not. lfw_as_salt_flx .and. liceform ) then SFWF_COMP(:,:,iblock,sfwf_comp_cpl) = FW(:,:,iblock) TFW_COMP (:,:,:,iblock,tfw_comp_cpl) = TFW(:,:,:,iblock) endif endif if ( luse_cpl_ifrac ) then OCN_WGT(:,:,iblock) = (c1-IFRAC(:,:,iblock)) * RCALCT(:,:,iblock) endif enddo !$OMP END PARALLEL DO !----------------------------------------------------------------------- ! Compute QSW_COSZ_WGHT_NORM. !----------------------------------------------------------------------- if ( qsw_distrb_iopt == qsw_distrb_iopt_cosz ) then tday00_interval_beg = tday00 !$OMP PARALLEL DO PRIVATE(iblock,nn,cosz_day) do iblock = 1, nblocks_clinic QSW_COSZ_WGHT_NORM(:,:,iblock) = c0 do nn = 1, nsteps_per_interval cosz_day = tday00_interval_beg + interval_cum_dayfrac(nn-1) & - interval_cum_dayfrac(nsteps_per_interval) call compute_cosz(cosz_day, iblock, QSW_COSZ_WGHT(:,:,iblock)) if (interval_avg_ts(nn)) then QSW_COSZ_WGHT_NORM(:,:,iblock) = & QSW_COSZ_WGHT_NORM(:,:,iblock) & + p5 * QSW_COSZ_WGHT(:,:,iblock) else QSW_COSZ_WGHT_NORM(:,:,iblock) = & QSW_COSZ_WGHT_NORM(:,:,iblock) & + QSW_COSZ_WGHT(:,:,iblock) endif enddo where (QSW_COSZ_WGHT_NORM(:,:,iblock) > c0) & QSW_COSZ_WGHT_NORM(:,:,iblock) = & (fullsteps_per_interval + p5 * halfsteps_per_interval) & / QSW_COSZ_WGHT_NORM(:,:,iblock) enddo !$OMP END PARALLEL DO endif #endif !----------------------------------------------------------------------- !EOC end subroutine pop_set_coupled_forcing !*********************************************************************** !BOP ! !IROUTINE: set_combined_forcing ! !INTERFACE: subroutine set_combined_forcing (STF,FW,TFW) ! !DESCRIPTION: ! ! This routine combines heat flux components into the STF array and ! converts from W/m**2, then combines terms when the "partially-coupled" ! has been selected ! ! !REVISION HISTORY: ! same as module ! !INPUT/OUTPUT PARAMETERS: real (r8), dimension(nx_block,ny_block,nt,max_blocks_clinic), & intent(inout) :: & STF, &! surface tracer fluxes at current timestep TFW ! tracer concentration in water flux real (r8), dimension(nx_block,ny_block,max_blocks_clinic), & intent(inout) :: & FW ! fresh water flux !EOP !BOC !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- integer (int_kind) :: & iblock, &! local address of current block n ! index #if CCSMCOUPLED real (r8), dimension(nx_block,ny_block,max_blocks_clinic) :: & WORK1, WORK2 ! local work arrays !*** need to zero out any padded cells WORK1 = c0 WORK2 = c0 if ( shf_formulation == 'partially-coupled' ) then !$OMP PARALLEL DO PRIVATE(iblock) do iblock=1,nblocks_clinic STF(:,:,1,iblock) = SHF_COMP(:,:,iblock,shf_comp_wrest) & + SHF_COMP(:,:,iblock,shf_comp_srest) & + SHF_COMP(:,:,iblock,shf_comp_cpl) enddo !$OMP END PARALLEL DO endif if ( sfwf_formulation == 'partially-coupled' ) then if (sfc_layer_type == sfc_layer_varthick .and. & .not. lfw_as_salt_flx) then !$OMP PARALLEL DO PRIVATE(iblock) do iblock=1,nblocks_clinic STF(:,:,2,iblock) = SFWF_COMP(:,:, iblock,sfwf_comp_wrest) & + SFWF_COMP(:,:, iblock,sfwf_comp_srest) FW(:,:,iblock) = SFWF_COMP(:,:, iblock,sfwf_comp_cpl) & + SFWF_COMP(:,:, iblock,sfwf_comp_flxio) TFW(:,:,:,iblock) = TFW_COMP(:,:,:,iblock, tfw_comp_cpl) & + TFW_COMP(:,:,:,iblock, tfw_comp_flxio) enddo !$OMP END PARALLEL DO else if ( lms_balance ) then !$OMP PARALLEL DO PRIVATE(iblock) do iblock=1,nblocks_clinic WORK1(:,:,iblock) = SFWF_COMP(:,:,iblock,sfwf_comp_flxio) / & salinity_factor WORK2(:,:,iblock) = SFWF_COMP(:,:,iblock,sfwf_comp_cpl) enddo !$OMP END PARALLEL DO call ms_balancing (WORK2, EVAP_F,PREC_F, MELT_F, ROFF_F, IOFF_F, & SALT_F, QFLUX, 'salt', ICEOCN_F=WORK1) !$OMP PARALLEL DO PRIVATE(iblock) do iblock=1,nblocks_clinic STF(:,:,2,iblock) = SFWF_COMP(:,:,iblock,sfwf_comp_wrest) & + SFWF_COMP(:,:,iblock,sfwf_comp_srest) & + WORK2(:,:,iblock) & + SFWF_COMP(:,:,iblock,sfwf_comp_flxio)* & MASK_SR(:,:,iblock) enddo !$OMP END PARALLEL DO else !$OMP PARALLEL DO PRIVATE(iblock) do iblock=1,nblocks_clinic STF(:,:,2,iblock) = SFWF_COMP(:,:,iblock,sfwf_comp_wrest) & + SFWF_COMP(:,:,iblock,sfwf_comp_srest) & + SFWF_COMP(:,:,iblock,sfwf_comp_cpl) & + SFWF_COMP(:,:,iblock,sfwf_comp_flxio) enddo !$OMP END PARALLEL DO endif endif endif #endif !----------------------------------------------------------------------- !EOC end subroutine set_combined_forcing !*********************************************************************** !BOP ! !IROUTINE: tavg_coupled_forcing ! !INTERFACE: subroutine tavg_coupled_forcing ! !DESCRIPTION: ! This routine accumulates tavg diagnostics related to forcing_coupled ! forcing. ! ! !REVISION HISTORY: ! same as module !EOP !BOC #if CCSMCOUPLED !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- integer (int_kind) :: & iblock ! block loop index type (block) :: & this_block ! block information for current block real (r8), dimension(nx_block,ny_block) :: & WORK ! local temp space for tavg diagnostics !----------------------------------------------------------------------- ! ! compute and accumulate tavg forcing diagnostics ! !----------------------------------------------------------------------- !$OMP PARALLEL DO PRIVATE(iblock,this_block) do iblock = 1,nblocks_clinic this_block = get_block(blocks_clinic(iblock),iblock) call accumulate_tavg_field(EVAP_F(:,:,iblock), tavg_EVAP_F,iblock,1) call accumulate_tavg_field(PREC_F(:,:,iblock), tavg_PREC_F,iblock,1) call accumulate_tavg_field(SNOW_F(:,:,iblock), tavg_SNOW_F,iblock,1) call accumulate_tavg_field(MELT_F(:,:,iblock), tavg_MELT_F,iblock,1) call accumulate_tavg_field(ROFF_F(:,:,iblock), tavg_ROFF_F,iblock,1) call accumulate_tavg_field(IOFF_F(:,:,iblock), tavg_IOFF_F,iblock,1) call accumulate_tavg_field(SALT_F(:,:,iblock), tavg_SALT_F,iblock,1) call accumulate_tavg_field(SENH_F(:,:,iblock), tavg_SENH_F,iblock,1) call accumulate_tavg_field(LWUP_F(:,:,iblock), tavg_LWUP_F,iblock,1) call accumulate_tavg_field(LWDN_F(:,:,iblock), tavg_LWDN_F,iblock,1) call accumulate_tavg_field(MELTH_F(:,:,iblock),tavg_MELTH_F,iblock,1) call accumulate_tavg_field(IFRAC(:,:,iblock), tavg_IFRAC,iblock,1) end do !$OMP END PARALLEL DO call tavg_mcog #endif !----------------------------------------------------------------------- !EOC end subroutine tavg_coupled_forcing !*********************************************************************** !BOP ! !IROUTINE: update_ghost_cells_coupler_fluxes ! !INTERFACE: subroutine update_ghost_cells_coupler_fluxes(errorCode) ! !DESCRIPTION: ! This routine accumulates tavg diagnostics related to forcing_coupled ! forcing. ! ! !REVISION HISTORY: ! same as module ! !OUTPUT PARAMETERS: integer (POP_i4), intent(out) :: errorCode !EOP !BOC !----------------------------------------------------------------------- ! ! update halos for all coupler fields ! !----------------------------------------------------------------------- errorCode = POP_Success #if CCSMCOUPLED call POP_HaloUpdate(SNOW_F,POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue = 0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'update_ghost_cells_coupler: error updating SNOW_F') return endif call POP_HaloUpdate(PREC_F,POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue = 0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'update_ghost_cells_coupler: error updating PREC_F') return endif call POP_HaloUpdate(EVAP_F,POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue = 0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'update_ghost_cells_coupler: error updating EVAP_F') return endif call POP_HaloUpdate(MELT_F,POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue = 0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'update_ghost_cells_coupler: error updating MELT_F') return endif call POP_HaloUpdate(ROFF_F,POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue = 0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'update_ghost_cells_coupler: error updating ROFF_F') return endif call POP_HaloUpdate(IOFF_F,POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue = 0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'update_ghost_cells_coupler: error updating IOFF_F') return endif call POP_HaloUpdate(SALT_F,POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue = 0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'update_ghost_cells_coupler: error updating SALT_F') return endif call POP_HaloUpdate(SENH_F,POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue = 0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'update_ghost_cells_coupler: error updating SENH_F') return endif call POP_HaloUpdate(LWUP_F,POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue = 0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'update_ghost_cells_coupler: error updating LWUP_F') return endif call POP_HaloUpdate(LWDN_F,POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue = 0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'update_ghost_cells_coupler: error updating LWDN_F') return endif call POP_HaloUpdate(MELTH_F,POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue = 0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'update_ghost_cells_coupler: error updating MELTH_F') return endif call POP_HaloUpdate(SHF_QSW,POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue = 0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'update_ghost_cells_coupler: error updating SHF_QSW') return endif call POP_HaloUpdate(IFRAC,POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue = 0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'update_ghost_cells_coupler: error updating IFRAC') return endif call POP_HaloUpdate(ATM_PRESS,POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue = 0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'update_ghost_cells_coupler: error updating ATM_PRESS') return endif call POP_HaloUpdate(U10_SQR,POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue = 0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'update_ghost_cells_coupler: error updating U10_SQR') return endif ! QL, 150526, LAMULT, USTOKES and VSTOKES call POP_HaloUpdate(LAMULT,POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue = 0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'update_ghost_cells_coupler: error updating LAMULT') return endif call POP_HaloUpdate(USTOKES,POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue = 0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'update_ghost_cells_coupler: error updating USTOKES') return endif call POP_HaloUpdate(VSTOKES,POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue = 0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'update_ghost_cells_coupler: error updating VSTOKES') return endif call POP_HaloUpdate(FRAC_BIN,POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue = 0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'update_ghost_cells_coupler: error updating FRAC_BIN') return endif call POP_HaloUpdate(QSW_RAW_BIN,POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue = 0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'update_ghost_cells_coupler: error updating QSW_RAW_BIN') return endif #endif !----------------------------------------------------------------------- !EOC end subroutine update_ghost_cells_coupler_fluxes !*********************************************************************** !BOP ! !IROUTINE: rotate_wind_stress ! !INTERFACE: subroutine rotate_wind_stress (WORK1,WORK2) ! !DESCRIPTION: ! This subroutine rotates true zonal/meridional wind stress into local ! coordinates, converts to dyne/cm**2, and shifts SMFT to the U grid ! ! !REVISION HISTORY: ! same as module ! !INPUT PARAMETERS: real (r8), dimension(nx_block,ny_block,max_blocks_clinic), intent(in) :: & WORK1, WORK2 ! contains taux and tauy from coupler !EOP !BOC #if CCSMCOUPLED !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- integer (kind=int_kind) :: iblock integer (POP_i4) :: errorCode !----------------------------------------------------------------------- ! ! rotate and convert ! !----------------------------------------------------------------------- SMFT(:,:,1,:) = (WORK1(:,:,:)*cos(ANGLET(:,:,:)) + & WORK2(:,:,:)*sin(ANGLET(:,:,:)))* & RCALCT(:,:,:)*momentum_factor SMFT(:,:,2,:) = (WORK2(:,:,:)*cos(ANGLET(:,:,:)) - & WORK1(:,:,:)*sin(ANGLET(:,:,:)))* & RCALCT(:,:,:)*momentum_factor !----------------------------------------------------------------------- ! ! perform halo updates following the vector rotation ! !----------------------------------------------------------------------- call POP_HaloUpdate(SMFT(:,:,1,:),POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindVector, errorCode, & fillValue = 0.0_POP_r8) call POP_HaloUpdate(SMFT(:,:,2,:),POP_haloClinic, & POP_gridHorzLocCenter, & POP_fieldKindVector, errorCode, & fillValue = 0.0_POP_r8) !----------------------------------------------------------------------- ! ! shift SMFT to U grid ! !----------------------------------------------------------------------- do iblock=1,nblocks_clinic call tgrid_to_ugrid(SMF(:,:,1,iblock),SMFT(:,:,1,iblock),iblock) call tgrid_to_ugrid(SMF(:,:,2,iblock),SMFT(:,:,2,iblock),iblock) enddo ! iblock #endif !----------------------------------------------------------------------- !EOC end subroutine rotate_wind_stress !*********************************************************************** !BOP ! !IROUTINE: compute_cosz ! !INTERFACE: subroutine compute_cosz(tday, iblock, COSZ) ! !DESCRIPTION: ! This subroutine computes cos of the solar zenith angle. ! Negative values are set to zero. ! ! !REVISION HISTORY: ! same as module ! ! !USES: use shr_orb_mod, only: shr_orb_decl, shr_orb_cosz ! !INPUT PARAMETERS: real (r8), intent(in) :: tday integer (int_kind), intent(in) :: iblock ! !OUTPUT PARAMETERS: real (r8), dimension(:,:), intent(out) :: COSZ !EOP !BOC !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- integer (int_kind) :: & i, j ! loop indices real (r8) :: & calday, & ! Calendar day, including fraction delta, & ! Solar declination angle in rad eccf ! Earth-sun distance factor (ie. (1/r)**2) !----------------------------------------------------------------------- call timer_start(timer_compute_cosz, block_id=iblock) ! shr_orb code assumes Jan 1 = calday 1, unlike Jan 1 = tday 0 calday = tday + c1 call shr_orb_decl(calday, orb_eccen, orb_mvelpp, orb_lambm0, & orb_obliqr, delta, eccf) do j = 1, ny_block do i = 1, nx_block COSZ(i,j) = shr_orb_cosz(calday, TLAT(i,j,iblock), & TLON(i,j,iblock), delta) COSZ(i,j) = max(c0, COSZ(i,j)) enddo enddo call timer_stop(timer_compute_cosz, block_id=iblock) !----------------------------------------------------------------------- !EOC end subroutine compute_cosz !*********************************************************************** end module forcing_coupled !|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||