! Copyright (c) 2004-2024 Lars Nerger ! ! This file is part of PDAF. ! ! PDAF is free software: you can redistribute it and/or modify ! it under the terms of the GNU Lesser General Public License ! as published by the Free Software Foundation, either version ! 3 of the License, or (at your option) any later version. ! ! PDAF is distributed in the hope that it will be useful, ! but WITHOUT ANY WARRANTY; without even the implied warranty of ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ! GNU Lesser General Public License for more details. ! ! You should have received a copy of the GNU Lesser General Public ! License along with PDAF. If not, see <http://www.gnu.org/licenses/>. ! !$Id$ !BOP ! ! !ROUTINE: PDAF_seek_update --- Control analysis update of the SEEK filter ! ! !INTERFACE: SUBROUTINE PDAF_seek_update(step, dim_p, dim_obs_p, dim_eof, state_p, & Uinv, eofV_p, eps, irediag, & U_init_dim_obs, U_obs_op, U_init_obs, U_prodRinvA, U_prepoststep, & screen, subtype, incremental, flag) ! !DESCRIPTION: ! Routine to control the analysis update of the SEEK filter. ! ! The analysis is performed by calling PDAF\_seek\_analysis ! and the rediagonalization is performed in PDAF\_seek\_rediag. ! In addition, the routine U\_prepoststep is called prior to ! the analysis and after the rediagonalization to allow the user ! to access the information of the modes and the state estimate. ! ! ! This is a core routine of PDAF and ! should not be changed by the user ! ! ! !REVISION HISTORY: ! 2003-07 - Lars Nerger - Initial code ! Later revisions - see svn log ! ! !USES: USE PDAF_timer, & ONLY: PDAF_timeit, PDAF_time_temp USE PDAF_mod_filtermpi, & ONLY: mype, dim_eof_l USE PDAF_mod_filter, & ONLY: forget, offline_mode IMPLICIT NONE ! !ARGUMENTS: INTEGER, INTENT(in) :: step ! Current time step INTEGER, INTENT(in) :: dim_p ! PE-local dimension of model state INTEGER, INTENT(out) :: dim_obs_p ! PE-local dimension of observation vector INTEGER, INTENT(in) :: dim_eof ! Number of EOFs REAL, INTENT(in) :: eps ! Epsilon for approximated TLM evolution INTEGER, INTENT(in) :: irediag ! Interval to perform rediagonalization REAL, INTENT(inout) :: state_p(dim_p) ! PE-local model state REAL, INTENT(inout) :: Uinv(dim_eof,dim_eof) ! Inverse of matrix U REAL, INTENT(inout) :: eofV_p(dim_p,dim_eof) ! PE-local matrix V INTEGER, INTENT(in) :: screen ! Verbosity flag INTEGER, INTENT(in) :: subtype ! Filter subtype INTEGER, INTENT(in) :: incremental ! Control incremental updating INTEGER, INTENT(inout) :: flag ! Status flag ! ! External subroutines ! ! (PDAF-internal names, real names are defined in the call to PDAF) EXTERNAL :: U_init_dim_obs, & ! Initialize dimension of observation vector U_obs_op, & ! Observation operator U_init_obs, & ! Initialize observation vector U_prepoststep, & ! User supplied pre/poststep routine U_prodRinvA ! Provide product R^-1 A ! !CALLING SEQUENCE: ! Called by: PDAF_put_state_seek ! Calls: U_prepoststep ! Calls: PDAF_seek_analysis ! Calls: PDAF_seek_rediag ! Calls: PDAF_timeit ! Calls: PDAF_time_temp !EOP ! *** local variables *** INTEGER :: i, col ! Counters INTEGER :: minusStep ! Time step counter REAL :: epsinv ! inverse of epsilon INTEGER :: countstep = 1 ! Internal step counter REAL, ALLOCATABLE :: Uinv_dyn(:,:) ! temporary matrix if Uinv is kept static ! ********************************************** ! *** Compute evolved basis of error space *** ! ********************************************** IF (subtype /= 2 .AND. subtype /= 3 .AND. .NOT.offline_mode) THEN ! Do not do mode-ensemble handling for fixed-basis variants epsinv = 1.0 / eps DO col = 1, dim_eof DO i = 1, dim_p eofV_p(i, col) = epsinv * (eofV_p(i, col) - state_p(i)) END DO END DO END IF ! ********************** ! *** Update phase *** ! ********************** ! *** prestep for forecast modes and state *** CALL PDAF_timeit(5, 'new') minusStep = -step ! Indicate forecast by negative time step number IF (mype == 0 .AND. screen > 0) THEN WRITE (*, '(a, 5x, a, i7)') 'PDAF', 'Call pre-post routine after forecast; step ', step ENDIF CALL U_prepoststep(minusStep, dim_p, dim_eof, dim_eof_l, dim_obs_p, & state_p, Uinv, eofV_p, flag) CALL PDAF_timeit(5, 'old') IF (mype == 0 .AND. screen > 0) THEN IF (screen > 1) THEN WRITE (*, '(a, 5x, a, F10.3, 1x, a)') & 'PDAF', '--- duration of prestep:', PDAF_time_temp(5), 's' END IF WRITE (*, '(a, 55a)') 'PDAF Analysis ', ('-', i = 1, 55) END IF #ifndef PDAF_NO_UPDATE CALL PDAF_timeit(3, 'new') ! *** SEEK analysis with forgetting factor *** subt: IF (subtype /= 3) THEN CALL PDAF_seek_analysis(step, dim_p, dim_obs_p, dim_eof, state_p, & Uinv, eofV_p, forget, U_init_dim_obs, U_obs_op, & U_init_obs, U_prodRinvA, screen, incremental, flag) ELSE subt ! For fixed Uinv initialize dynamic Uinv which is hold only here ALLOCATE(Uinv_dyn(dim_eof, dim_eof)) Uinv_dyn = Uinv ! Perform analysis CALL PDAF_seek_analysis(step, dim_p, dim_obs_p, dim_eof, state_p, & Uinv_dyn, eofV_p, forget, U_init_dim_obs, U_obs_op, & U_init_obs, U_prodRinvA, screen, incremental, flag) END IF subt CALL PDAF_timeit(3, 'old') IF (mype == 0 .AND. screen > 1) THEN WRITE (*, '(a, 5x, a, F10.3, 1x, a)') & 'PDAF', '--- update duration:', PDAF_time_temp(3), 's' END IF ! *** Re-orthogonalize the covariance modes *** re_diag: IF (irediag > 0) THEN re_diag2: IF ( MOD(countstep, irediag) == 0) THEN CALL PDAF_timeit(4, 'new') CALL PDAF_seek_rediag(dim_p, dim_eof, Uinv, eofV_p, subtype, & screen, flag) CALL PDAF_timeit(4, 'old') IF (mype == 0 .AND. screen > 1) THEN WRITE (*, '(a, 5x, a, F10.3, 1x, a)') & 'PDAF', '--- re-diag duration:', PDAF_time_temp(4), 's' END IF END IF re_diag2 END IF re_diag ! increment stepping counter countstep = countstep + 1 #else WRITE (*,'(/5x,a/)') & '!!! PDAF WARNING: ANALYSIS STEP IS DEACTIVATED BY PDAF_NO_UPDATE !!!' #endif ! *** poststep for analysis modes and state *** CALL PDAF_timeit(5, 'new') IF (mype == 0 .AND. screen > 0) THEN WRITE (*, '(a, 5x, a)') 'PDAF', 'Call pre-post routine after analysis step' ENDIF IF (subtype /= 3) THEN CALL U_prepoststep(step, dim_p, dim_eof, dim_eof_l, dim_obs_p, & state_p, Uinv, eofV_p, flag) ELSE ! prepoststep with fixed Uinv - hand over Uinv as changed by analysis CALL U_prepoststep(step, dim_p, dim_eof, dim_eof_l, dim_obs_p, & state_p, Uinv_dyn, eofV_p, flag) DEALLOCATE(Uinv_dyn) END IF CALL PDAF_timeit(5, 'old') IF (mype == 0 .AND. screen > 0) THEN IF (screen > 1) THEN WRITE (*, '(a, 5x, a, F10.3, 1x, a)') & 'PDAF', '--- duration of poststep:', PDAF_time_temp(5), 's' END IF WRITE (*, '(a, 55a)') 'PDAF Forecast ', ('-', i = 1, 55) END IF END SUBROUTINE PDAF_seek_update