Changeset 11332


Ignore:
Timestamp:
2019-07-23T16:26:43+02:00 (13 months ago)
Author:
clem
Message:

prepare code to read temperature and salinity for bdy

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icevar.F90

    r11223 r11332  
    773773   !! ** Purpose :  converting N-cat ice to jpl ice categories 
    774774   !!------------------------------------------------------------------- 
    775    SUBROUTINE ice_var_itd_1c1c( zhti, zhts, zati, zh_i, zh_s, za_i ) 
     775   SUBROUTINE ice_var_itd_1c1c( phti, phts, pati ,       ph_i, ph_s, pa_i, & 
     776      &                         ptmi, ptms, ptmsu, psmi, pt_i, pt_s, pt_su, ps_i ) 
    776777      !!------------------------------------------------------------------- 
    777778      !! ** Purpose :  converting 1-cat ice to 1 ice category 
    778779      !!------------------------------------------------------------------- 
    779       REAL(wp), DIMENSION(:), INTENT(in)    ::   zhti, zhts, zati    ! input ice/snow variables 
    780       REAL(wp), DIMENSION(:), INTENT(inout) ::   zh_i, zh_s, za_i    ! output ice/snow variables 
    781       !!------------------------------------------------------------------- 
    782       zh_i(:) = zhti(:) 
    783       zh_s(:) = zhts(:) 
    784       za_i(:) = zati(:) 
     780      REAL(wp), DIMENSION(:), INTENT(in)    ::   phti, phts, pati    ! input  ice/snow variables 
     781      REAL(wp), DIMENSION(:), INTENT(inout) ::   ph_i, ph_s, pa_i    ! output ice/snow variables 
     782      REAL(wp), DIMENSION(:), INTENT(in)   , OPTIONAL ::   ptmi, ptms, ptmsu, psmi    ! input  ice/snow temp & sal 
     783      REAL(wp), DIMENSION(:), INTENT(inout), OPTIONAL ::   pt_i, pt_s, pt_su, ps_i    ! output ice/snow temp & sal 
     784      !!------------------------------------------------------------------- 
     785      ! == thickness and concentration == ! 
     786      ph_i(:) = phti(:) 
     787      ph_s(:) = phts(:) 
     788      pa_i(:) = pati(:) 
     789      ! 
     790      ! == temperature and salinity == ! 
     791      IF( PRESENT( pt_i  ) )   pt_i (:) = ptmi (:) 
     792      IF( PRESENT( pt_s  ) )   pt_s (:) = ptms (:) 
     793      IF( PRESENT( pt_su ) )   pt_su(:) = ptmsu(:) 
     794      IF( PRESENT( ps_i  ) )   ps_i (:) = psmi (:) 
     795       
    785796   END SUBROUTINE ice_var_itd_1c1c 
    786797 
    787    SUBROUTINE ice_var_itd_Nc1c( zhti, zhts, zati, zh_i, zh_s, za_i ) 
     798   SUBROUTINE ice_var_itd_Nc1c( phti, phts, pati ,       ph_i, ph_s, pa_i, & 
     799      &                         ptmi, ptms, ptmsu, psmi, pt_i, pt_s, pt_su, ps_i ) 
    788800      !!------------------------------------------------------------------- 
    789801      !! ** Purpose :  converting N-cat ice to 1 ice category 
    790802      !!------------------------------------------------------------------- 
    791       REAL(wp), DIMENSION(:,:), INTENT(in)    ::   zhti, zhts, zati    ! input ice/snow variables 
    792       REAL(wp), DIMENSION(:)  , INTENT(inout) ::   zh_i, zh_s, za_i    ! output ice/snow variables 
    793       !!------------------------------------------------------------------- 
    794       ! 
    795       za_i(:) = SUM( zati(:,:), dim=2 ) 
    796       ! 
    797       WHERE( za_i(:) /= 0._wp ) 
    798          zh_i(:) = SUM( zhti(:,:) * zati(:,:), dim=2 ) / za_i(:) 
    799          zh_s(:) = SUM( zhts(:,:) * zati(:,:), dim=2 ) / za_i(:) 
    800       ELSEWHERE 
    801          zh_i(:) = 0._wp 
    802          zh_s(:) = 0._wp 
     803      REAL(wp), DIMENSION(:,:), INTENT(in)    ::   phti, phts, pati    ! input  ice/snow variables 
     804      REAL(wp), DIMENSION(:)  , INTENT(inout) ::   ph_i, ph_s, pa_i    ! output ice/snow variables 
     805      REAL(wp), DIMENSION(:,:), INTENT(in)   , OPTIONAL ::   ptmi, ptms, ptmsu, psmi    ! input  ice/snow temp & sal 
     806      REAL(wp), DIMENSION(:)  , INTENT(inout), OPTIONAL ::   pt_i, pt_s, pt_su, ps_i    ! output ice/snow temp & sal 
     807      ! 
     808      REAL(wp), ALLOCATABLE, DIMENSION(:) ::   z1_ai, z1_vi, z1_vs 
     809      ! 
     810      INTEGER ::   idim   
     811      !!------------------------------------------------------------------- 
     812      ! 
     813      idim = SIZE( phti, 1 ) 
     814      ! 
     815      ! == thickness and concentration == ! 
     816      ALLOCATE( z1_ai(idim) ) 
     817      ! 
     818      pa_i(:) = SUM( pati(:,:), dim=2 ) 
     819 
     820      WHERE( ( pa_i(:) ) /= 0._wp )   ;   z1_ai(:) = 1._wp / pa_i(:) 
     821      ELSEWHERE                       ;   z1_ai(:) = 0._wp 
    803822      END WHERE 
     823 
     824      ph_i(:) = SUM( phti(:,:) * pati(:,:), dim=2 ) * z1_ai(:) 
     825      ph_s(:) = SUM( phts(:,:) * pati(:,:), dim=2 ) * z1_ai(:) 
     826      ! 
     827      ! == temperature and salinity == ! 
     828      IF( PRESENT( pt_i ) .OR. PRESENT( pt_s ) .OR. PRESENT( pt_su ) .OR. PRESENT( ps_i ) ) THEN 
     829         ! 
     830         ALLOCATE( z1_vi(idim), z1_vs(idim) ) 
     831         ! 
     832         WHERE( ( pa_i(:) * ph_i(:) ) /= 0._wp )   ;   z1_vi(:) = 1._wp / ( pa_i(:) * ph_i(:) ) 
     833         ELSEWHERE                                 ;   z1_vi(:) = 0._wp 
     834         END WHERE 
     835         WHERE( ( pa_i(:) * ph_s(:) ) /= 0._wp )   ;   z1_vs(:) = 1._wp / ( pa_i(:) * ph_s(:) ) 
     836         ELSEWHERE                                 ;   z1_vs(:) = 0._wp 
     837         END WHERE 
     838         ! 
     839         IF( PRESENT( pt_i  ) )   pt_i (:) = SUM( ptmi (:,:) * pati(:,:) * phti(:,:), dim=2 ) * z1_vi(:) 
     840         IF( PRESENT( pt_s  ) )   pt_s (:) = SUM( ptms (:,:) * pati(:,:) * phts(:,:), dim=2 ) * z1_vs(:) 
     841         IF( PRESENT( pt_su ) )   pt_su(:) = SUM( ptmsu(:,:) * pati(:,:)            , dim=2 ) * z1_ai(:) 
     842         IF( PRESENT( ps_i  ) )   ps_i (:) = SUM( psmi (:,:) * pati(:,:) * phti(:,:), dim=2 ) * z1_vi(:) 
     843         ! 
     844         DEALLOCATE( z1_vi, z1_vs ) 
     845         ! 
     846      ENDIF 
     847      ! 
     848      DEALLOCATE( z1_ai ) 
    804849      ! 
    805850   END SUBROUTINE ice_var_itd_Nc1c 
    806851    
    807    SUBROUTINE ice_var_itd_1cMc( zhti, zhts, zati, zh_i, zh_s, za_i ) 
     852   SUBROUTINE ice_var_itd_1cMc( phti, phts, pati ,       ph_i, ph_s, pa_i, & 
     853      &                         ptmi, ptms, ptmsu, psmi, pt_i, pt_s, pt_su, ps_i ) 
    808854      !!------------------------------------------------------------------- 
    809855      !! 
     
    826872      !!               4) Iterate until ok (SUM(itest(:) = 4) 
    827873      !! 
    828       !! ** Arguments : zhti: 1-cat ice thickness 
    829       !!                zhts: 1-cat snow depth 
    830       !!                zati: 1-cat ice concentration 
     874      !! ** Arguments : phti: 1-cat ice thickness 
     875      !!                phts: 1-cat snow depth 
     876      !!                pati: 1-cat ice concentration 
    831877      !! 
    832878      !! ** Output    : jpl-cat  
     
    834880      !!  (Example of application: BDY forcings when input are cell averaged)   
    835881      !!------------------------------------------------------------------- 
    836       INTEGER  :: ji, jk, jl             ! dummy loop indices 
    837       INTEGER  :: idim, i_fill, jl0   
    838       REAL(wp) :: zarg, zV, zconv, zdh, zdv 
    839       REAL(wp), DIMENSION(:),   INTENT(in)    ::   zhti, zhts, zati    ! input  ice/snow variables 
    840       REAL(wp), DIMENSION(:,:), INTENT(inout) ::   zh_i, zh_s, za_i    ! output ice/snow variables 
    841       INTEGER , DIMENSION(4)                  ::   itest 
    842       !!------------------------------------------------------------------- 
    843       ! 
    844       ! ---------------------------------------- 
     882      REAL(wp), DIMENSION(:),   INTENT(in)    ::   phti, phts, pati    ! input  ice/snow variables 
     883      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   ph_i, ph_s, pa_i    ! output ice/snow variables 
     884      REAL(wp), DIMENSION(:)  , INTENT(in)   , OPTIONAL ::   ptmi, ptms, ptmsu, psmi    ! input  ice/snow temp & sal 
     885      REAL(wp), DIMENSION(:,:), INTENT(inout), OPTIONAL ::   pt_i, pt_s, pt_su, ps_i    ! output ice/snow temp & sal 
     886      ! 
     887      INTEGER , DIMENSION(4) ::   itest 
     888      INTEGER  ::   ji, jk, jl 
     889      INTEGER  ::   idim, i_fill, jl0   
     890      REAL(wp) ::   zarg, zV, zconv, zdh, zdv 
     891      !!------------------------------------------------------------------- 
     892      ! 
     893      ! == thickness and concentration == ! 
    845894      ! distribution over the jpl ice categories 
    846       ! ---------------------------------------- 
    847       ! a gaussian distribution for ice concentration is used 
    848       ! then we check whether the distribution fullfills 
    849       ! volume and area conservation, positivity and ice categories bounds 
    850       idim = SIZE( zhti , 1 ) 
    851       zh_i(1:idim,1:jpl) = 0._wp 
    852       zh_s(1:idim,1:jpl) = 0._wp 
    853       za_i(1:idim,1:jpl) = 0._wp 
    854       ! 
    855       IF( jpl == 1 ) THEN 
    856          CALL ice_var_itd_1c1c( zhti, zhts, zati, zh_i(:,1), zh_s(:,1), za_i(:,1) ) 
    857          RETURN 
    858       ENDIF 
     895      !    a gaussian distribution for ice concentration is used 
     896      !    then we check whether the distribution fullfills 
     897      !    volume and area conservation, positivity and ice categories bounds 
     898      idim = SIZE( phti , 1 ) 
     899      ! 
     900      ph_i(1:idim,1:jpl) = 0._wp 
     901      ph_s(1:idim,1:jpl) = 0._wp 
     902      pa_i(1:idim,1:jpl) = 0._wp 
    859903      ! 
    860904      DO ji = 1, idim 
    861905         ! 
    862          IF( zhti(ji) > 0._wp ) THEN 
     906         IF( phti(ji) > 0._wp ) THEN 
    863907            ! 
    864908            ! find which category (jl0) the input ice thickness falls into 
    865909            jl0 = jpl 
    866910            DO jl = 1, jpl 
    867                IF ( ( zhti(ji) >= hi_max(jl-1) ) .AND. ( zhti(ji) < hi_max(jl) ) ) THEN 
     911               IF ( ( phti(ji) >= hi_max(jl-1) ) .AND. ( phti(ji) < hi_max(jl) ) ) THEN 
    868912                  jl0 = jl 
    869913                  CYCLE 
     
    877921               i_fill = i_fill - 1 
    878922               ! 
    879                zh_i(ji,1:jpl) = 0._wp 
    880                za_i(ji,1:jpl) = 0._wp 
     923               ph_i(ji,1:jpl) = 0._wp 
     924               pa_i(ji,1:jpl) = 0._wp 
    881925               itest(:)       = 0       
    882926               ! 
    883927               IF ( i_fill == 1 ) THEN      !-- case very thin ice: fill only category 1 
    884                   zh_i(ji,1) = zhti(ji) 
    885                   za_i (ji,1) = zati (ji) 
     928                  ph_i(ji,1) = phti(ji) 
     929                  pa_i (ji,1) = pati (ji) 
    886930               ELSE                         !-- case ice is thicker: fill categories >1 
    887931                  ! thickness 
    888932                  DO jl = 1, i_fill - 1 
    889                      zh_i(ji,jl) = hi_mean(jl) 
     933                     ph_i(ji,jl) = hi_mean(jl) 
    890934                  END DO 
    891935                  ! 
    892936                  ! concentration 
    893                   za_i(ji,jl0) = zati(ji) / SQRT(REAL(jpl)) 
     937                  pa_i(ji,jl0) = pati(ji) / SQRT(REAL(jpl)) 
    894938                  DO jl = 1, i_fill - 1 
    895939                     IF ( jl /= jl0 ) THEN 
    896                         zarg        = ( zh_i(ji,jl) - zhti(ji) ) / ( zhti(ji) * 0.5_wp ) 
    897                         za_i(ji,jl) =   za_i (ji,jl0) * EXP(-zarg**2) 
     940                        zarg        = ( ph_i(ji,jl) - phti(ji) ) / ( phti(ji) * 0.5_wp ) 
     941                        pa_i(ji,jl) =   pa_i (ji,jl0) * EXP(-zarg**2) 
    898942                     ENDIF 
    899943                  END DO 
    900944                  ! 
    901945                  ! last category 
    902                   za_i(ji,i_fill) = zati(ji) - SUM( za_i(ji,1:i_fill-1) ) 
    903                   zV = SUM( za_i(ji,1:i_fill-1) * zh_i(ji,1:i_fill-1) ) 
    904                   zh_i(ji,i_fill) = ( zhti(ji) * zati(ji) - zV ) / MAX( za_i(ji,i_fill), epsi10 )  
     946                  pa_i(ji,i_fill) = pati(ji) - SUM( pa_i(ji,1:i_fill-1) ) 
     947                  zV = SUM( pa_i(ji,1:i_fill-1) * ph_i(ji,1:i_fill-1) ) 
     948                  ph_i(ji,i_fill) = ( phti(ji) * pati(ji) - zV ) / MAX( pa_i(ji,i_fill), epsi10 )  
    905949                  ! 
    906950                  ! correction if concentration of upper cat is greater than lower cat 
     
    908952                  IF ( jl0 /= jpl ) THEN 
    909953                     DO jl = jpl, jl0+1, -1 
    910                         IF ( za_i(ji,jl) > za_i(ji,jl-1) ) THEN 
    911                            zdv = zh_i(ji,jl) * za_i(ji,jl) 
    912                            zh_i(ji,jl    ) = 0._wp 
    913                            za_i (ji,jl    ) = 0._wp 
    914                            za_i (ji,1:jl-1) = za_i(ji,1:jl-1) + zdv / MAX( REAL(jl-1) * zhti(ji), epsi10 ) 
     954                        IF ( pa_i(ji,jl) > pa_i(ji,jl-1) ) THEN 
     955                           zdv = ph_i(ji,jl) * pa_i(ji,jl) 
     956                           ph_i(ji,jl    ) = 0._wp 
     957                           pa_i (ji,jl    ) = 0._wp 
     958                           pa_i (ji,1:jl-1) = pa_i(ji,1:jl-1) + zdv / MAX( REAL(jl-1) * phti(ji), epsi10 ) 
    915959                        END IF 
    916960                     END DO 
     
    920964               ! 
    921965               ! Compatibility tests 
    922                zconv = ABS( zati(ji) - SUM( za_i(ji,1:jpl) ) )  
     966               zconv = ABS( pati(ji) - SUM( pa_i(ji,1:jpl) ) )  
    923967               IF ( zconv < epsi06 )   itest(1) = 1                                        ! Test 1: area conservation 
    924968               ! 
    925                zconv = ABS( zhti(ji)*zati(ji) - SUM( za_i(ji,1:jpl)*zh_i(ji,1:jpl) ) ) 
     969               zconv = ABS( phti(ji)*pati(ji) - SUM( pa_i(ji,1:jpl)*ph_i(ji,1:jpl) ) ) 
    926970               IF ( zconv < epsi06 )   itest(2) = 1                                        ! Test 2: volume conservation 
    927971               ! 
    928                IF ( zh_i(ji,i_fill) >= hi_max(i_fill-1) )   itest(3) = 1                  ! Test 3: thickness of the last category is in-bounds ? 
     972               IF ( ph_i(ji,i_fill) >= hi_max(i_fill-1) )   itest(3) = 1                  ! Test 3: thickness of the last category is in-bounds ? 
    929973               ! 
    930974               itest(4) = 1 
    931975               DO jl = 1, i_fill 
    932                   IF ( za_i(ji,jl) < 0._wp ) itest(4) = 0                                ! Test 4: positivity of ice concentrations 
     976                  IF ( pa_i(ji,jl) < 0._wp ) itest(4) = 0                                ! Test 4: positivity of ice concentrations 
    933977               END DO 
    934978               !                                         !---------------------------- 
     
    938982      END DO 
    939983 
    940       ! Add Snow in each category where za_i is not 0 
     984      ! Add Snow in each category where pa_i is not 0 
    941985      DO jl = 1, jpl 
    942986         DO ji = 1, idim 
    943             IF( za_i(ji,jl) > 0._wp ) THEN 
    944                zh_s(ji,jl) = zh_i(ji,jl) * ( zhts(ji) / zhti(ji) ) 
     987            IF( pa_i(ji,jl) > 0._wp ) THEN 
     988               ph_s(ji,jl) = ph_i(ji,jl) * ( phts(ji) / phti(ji) ) 
    945989               ! In case snow load is in excess that would lead to transformation from snow to ice 
    946990               ! Then, transfer the snow excess into the ice (different from icethd_dh) 
    947                zdh = MAX( 0._wp, ( rhos * zh_s(ji,jl) + ( rhoi - rau0 ) * zh_i(ji,jl) ) * r1_rau0 )  
     991               zdh = MAX( 0._wp, ( rhos * ph_s(ji,jl) + ( rhoi - rau0 ) * ph_i(ji,jl) ) * r1_rau0 )  
    948992               ! recompute h_i, h_s avoiding out of bounds values 
    949                zh_i(ji,jl) = MIN( hi_max(jl), zh_i(ji,jl) + zdh ) 
    950                zh_s(ji,jl) = MAX( 0._wp, zh_s(ji,jl) - zdh * rhoi * r1_rhos ) 
     993               ph_i(ji,jl) = MIN( hi_max(jl), ph_i(ji,jl) + zdh ) 
     994               ph_s(ji,jl) = MAX( 0._wp, ph_s(ji,jl) - zdh * rhoi * r1_rhos ) 
    951995            ENDIF 
    952996         END DO 
    953997      END DO 
    954998      ! 
     999      ! == temperature and salinity == ! 
     1000      IF( PRESENT( pt_i  ) ) THEN 
     1001         DO jl = 1, jpl 
     1002            pt_i(:,jl) = ptmi (:) 
     1003         END DO 
     1004      ENDIF 
     1005      IF( PRESENT( pt_s  ) ) THEN 
     1006         DO jl = 1, jpl 
     1007            pt_s (:,jl) = ptms (:) 
     1008         END DO 
     1009      ENDIF 
     1010      IF( PRESENT( pt_su ) ) THEN 
     1011         DO jl = 1, jpl 
     1012            pt_su(:,jl) = ptmsu(:) 
     1013         END DO 
     1014      ENDIF 
     1015      IF( PRESENT( ps_i  ) ) THEN 
     1016         DO jl = 1, jpl 
     1017            ps_i (:,jl) = psmi (:) 
     1018         END DO 
     1019      ENDIF 
     1020      ! 
    9551021   END SUBROUTINE ice_var_itd_1cMc 
    9561022 
    957    SUBROUTINE ice_var_itd_NcMc( zhti, zhts, zati, zh_i, zh_s, za_i ) 
     1023   SUBROUTINE ice_var_itd_NcMc( phti, phts, pati ,       ph_i, ph_s, pa_i, & 
     1024      &                         ptmi, ptms, ptmsu, psmi, pt_i, pt_s, pt_su, ps_i ) 
    9581025      !!------------------------------------------------------------------- 
    9591026      !! 
     
    9761043      !!                      b) removing 25% ice area from the higher cat (descendant loop jlmax=>jlmin) 
    9771044      !! 
    978       !! ** Arguments : zhti: N-cat ice thickness 
    979       !!                zhts: N-cat snow depth 
    980       !!                zati: N-cat ice concentration 
     1045      !! ** Arguments : phti: N-cat ice thickness 
     1046      !!                phts: N-cat snow depth 
     1047      !!                pati: N-cat ice concentration 
    9811048      !! 
    9821049      !! ** Output    : jpl-cat  
     
    9841051      !!  (Example of application: BDY forcings when inputs have N-cat /= jpl)   
    9851052      !!------------------------------------------------------------------- 
    986       INTEGER  ::   ji, jl, jl1, jl2             ! dummy loop indices 
     1053      REAL(wp), DIMENSION(:,:), INTENT(in)    ::   phti, phts, pati    ! input  ice/snow variables 
     1054      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   ph_i, ph_s, pa_i    ! output ice/snow variables 
     1055      REAL(wp), DIMENSION(:,:), INTENT(in)   , OPTIONAL ::   ptmi, ptms, ptmsu, psmi    ! input  ice/snow temp & sal 
     1056      REAL(wp), DIMENSION(:,:), INTENT(inout), OPTIONAL ::   pt_i, pt_s, pt_su, ps_i    ! output ice/snow temp & sal 
     1057      ! 
     1058      INTEGER , ALLOCATABLE, DIMENSION(:,:) ::   jlfil, jlfil2 
     1059      INTEGER , ALLOCATABLE, DIMENSION(:)   ::   jlmax, jlmin 
     1060      REAL(wp), ALLOCATABLE, DIMENSION(:)   ::   z1_ai, z1_vi, z1_vs, ztmp 
     1061      ! 
     1062      REAL(wp), PARAMETER ::   ztrans = 0.25_wp 
     1063      INTEGER  ::   ji, jl, jl1, jl2 
    9871064      INTEGER  ::   idim, icat   
    988       REAL(wp), PARAMETER ::   ztrans = 0.25_wp 
    989       REAL(wp), DIMENSION(:,:), INTENT(in)    ::   zhti, zhts, zati    ! input ice/snow variables 
    990       REAL(wp), DIMENSION(:,:), INTENT(inout) ::   zh_i, zh_s, za_i    ! output ice/snow variables 
    991       INTEGER , DIMENSION(:,:), ALLOCATABLE   ::   jlfil, jlfil2 
    992       INTEGER , DIMENSION(:)  , ALLOCATABLE   ::   jlmax, jlmin 
    993       !!------------------------------------------------------------------- 
    994       ! 
    995       idim = SIZE( zhti, 1 ) 
    996       icat = SIZE( zhti, 2 ) 
     1065      !!------------------------------------------------------------------- 
     1066      ! 
     1067      idim = SIZE( phti, 1 ) 
     1068      icat = SIZE( phti, 2 ) 
     1069      ! 
     1070      ! == thickness and concentration == ! 
    9971071      !                                 ! ---------------------- ! 
    9981072      IF( icat == jpl ) THEN            ! input cat = output cat ! 
    9991073         !                              ! ---------------------- ! 
    1000          zh_i(:,:) = zhti(:,:) 
    1001          zh_s(:,:) = zhts(:,:) 
    1002          za_i(:,:) = zati(:,:) 
     1074         ph_i(:,:) = phti(:,:) 
     1075         ph_s(:,:) = phts(:,:) 
     1076         pa_i(:,:) = pati(:,:) 
     1077         ! 
     1078         ! == temperature and salinity == ! 
     1079         IF( PRESENT( pt_i  ) )   pt_i (:,:) = ptmi (:,:) 
     1080         IF( PRESENT( pt_s  ) )   pt_s (:,:) = ptms (:,:) 
     1081         IF( PRESENT( pt_su ) )   pt_su(:,:) = ptmsu(:,:) 
     1082         IF( PRESENT( ps_i  ) )   ps_i (:,:) = psmi (:,:) 
    10031083         !                              ! ---------------------- ! 
    1004       ELSEIF( icat == 1 ) THEN          ! specific case if N = 1 ! 
     1084      ELSEIF( icat == 1 ) THEN          ! input cat = 1          ! 
    10051085         !                              ! ---------------------- ! 
    1006          ! 
    1007          CALL ice_var_itd_1cMc( zhti(:,1), zhts(:,1), zati(:,1), zh_i, zh_s, za_i ) 
    1008          ! 
     1086         CALL  ice_var_itd_1cMc( phti(:,1), phts(:,1), pati (:,1),            ph_i(:,:), ph_s(:,:), pa_i (:,:), & 
     1087            &                    ptmi(:,1), ptms(:,1), ptmsu(:,1), psmi(:,1), pt_i(:,:), pt_s(:,:), pt_su(:,:), ps_i(:,:) ) 
    10091088         !                              ! ---------------------- ! 
    1010       ELSEIF( jpl == 1 ) THEN           ! specific case if M = 1 ! 
     1089      ELSEIF( jpl == 1 ) THEN           ! output cat = 1        ! 
    10111090         !                              ! ---------------------- ! 
    1012          ! 
    1013          CALL ice_var_itd_Nc1c( zhti, zhts, zati, zh_i(:,1), zh_s(:,1), za_i(:,1) ) 
    1014          ! 
     1091         CALL  ice_var_itd_Nc1c( phti(:,:), phts(:,:), pati (:,:),            ph_i(:,1), ph_s(:,1), pa_i (:,1), & 
     1092            &                    ptmi(:,:), ptms(:,:), ptmsu(:,:), psmi(:,:), pt_i(:,1), pt_s(:,1), pt_su(:,1), ps_i(:,1) )          
    10151093         !                              ! ----------------------- ! 
    10161094      ELSE                              ! input cat /= output cat ! 
     
    10211099 
    10221100         ! --- initialize output fields to 0 --- ! 
    1023          zh_i(1:idim,1:jpl) = 0._wp 
    1024          zh_s(1:idim,1:jpl) = 0._wp 
    1025          za_i(1:idim,1:jpl) = 0._wp 
     1101         ph_i(1:idim,1:jpl) = 0._wp 
     1102         ph_s(1:idim,1:jpl) = 0._wp 
     1103         pa_i(1:idim,1:jpl) = 0._wp 
    10261104         ! 
    10271105         ! --- fill the categories --- ! 
     
    10331111            DO jl2 = 1, icat 
    10341112               DO ji = 1, idim 
    1035                   IF( hi_max(jl1-1) <= zhti(ji,jl2) .AND. hi_max(jl1) > zhti(ji,jl2) ) THEN 
     1113                  IF( hi_max(jl1-1) <= phti(ji,jl2) .AND. hi_max(jl1) > phti(ji,jl2) ) THEN 
    10361114                     ! fill the right category 
    1037                      zh_i(ji,jl1) = zhti(ji,jl2) 
    1038                      zh_s(ji,jl1) = zhts(ji,jl2) 
    1039                      za_i(ji,jl1) = zati(ji,jl2) 
     1115                     ph_i(ji,jl1) = phti(ji,jl2) 
     1116                     ph_s(ji,jl1) = phts(ji,jl2) 
     1117                     pa_i(ji,jl1) = pati(ji,jl2) 
    10401118                     ! record categories that are filled 
    10411119                     jlmax(ji) = MAX( jlmax(ji), jl1 ) 
     
    10541132            IF( jl1 > 1 ) THEN 
    10551133               ! fill the lower cat (jl1-1) 
    1056                za_i(ji,jl1-1) = ztrans * za_i(ji,jl1) 
    1057                zh_i(ji,jl1-1) = hi_mean(jl1-1) 
     1134               pa_i(ji,jl1-1) = ztrans * pa_i(ji,jl1) 
     1135               ph_i(ji,jl1-1) = hi_mean(jl1-1) 
    10581136               ! remove from cat jl1 
    1059                za_i(ji,jl1  ) = ( 1._wp - ztrans ) * za_i(ji,jl1) 
     1137               pa_i(ji,jl1  ) = ( 1._wp - ztrans ) * pa_i(ji,jl1) 
    10601138            ENDIF 
    10611139            IF( jl2 < jpl ) THEN 
    10621140               ! fill the upper cat (jl2+1) 
    1063                za_i(ji,jl2+1) = ztrans * za_i(ji,jl2) 
    1064                zh_i(ji,jl2+1) = hi_mean(jl2+1) 
     1141               pa_i(ji,jl2+1) = ztrans * pa_i(ji,jl2) 
     1142               ph_i(ji,jl2+1) = hi_mean(jl2+1) 
    10651143               ! remove from cat jl2 
    1066                za_i(ji,jl2  ) = ( 1._wp - ztrans ) * za_i(ji,jl2) 
     1144               pa_i(ji,jl2  ) = ( 1._wp - ztrans ) * pa_i(ji,jl2) 
    10671145            ENDIF 
    10681146         END DO 
     
    10741152               IF( jlfil(ji,jl-1) /= 0 .AND. jlfil(ji,jl) == 0 ) THEN 
    10751153                  ! fill high 
    1076                   za_i(ji,jl) = ztrans * za_i(ji,jl-1) 
    1077                   zh_i(ji,jl) = hi_mean(jl) 
     1154                  pa_i(ji,jl) = ztrans * pa_i(ji,jl-1) 
     1155                  ph_i(ji,jl) = hi_mean(jl) 
    10781156                  jlfil(ji,jl) = jl 
    10791157                  ! remove low 
    1080                   za_i(ji,jl-1) = ( 1._wp - ztrans ) * za_i(ji,jl-1) 
     1158                  pa_i(ji,jl-1) = ( 1._wp - ztrans ) * pa_i(ji,jl-1) 
    10811159               ENDIF 
    10821160            END DO 
     
    10881166               IF( jlfil2(ji,jl+1) /= 0 .AND. jlfil2(ji,jl) == 0 ) THEN 
    10891167                  ! fill low 
    1090                   za_i(ji,jl) = za_i(ji,jl) + ztrans * za_i(ji,jl+1) 
    1091                   zh_i(ji,jl) = hi_mean(jl)  
     1168                  pa_i(ji,jl) = pa_i(ji,jl) + ztrans * pa_i(ji,jl+1) 
     1169                  ph_i(ji,jl) = hi_mean(jl)  
    10921170                  jlfil2(ji,jl) = jl 
    10931171                  ! remove high 
    1094                   za_i(ji,jl+1) = ( 1._wp - ztrans ) * za_i(ji,jl+1) 
     1172                  pa_i(ji,jl+1) = ( 1._wp - ztrans ) * pa_i(ji,jl+1) 
    10951173               ENDIF 
    10961174            END DO 
     
    10991177         DEALLOCATE( jlfil, jlfil2 )      ! deallocate arrays 
    11001178         DEALLOCATE( jlmin, jlmax ) 
     1179         ! 
     1180         ! == temperature and salinity == ! 
     1181         ! 
     1182         IF( PRESENT( pt_i ) .OR. PRESENT( pt_s ) .OR. PRESENT( pt_su ) .OR. PRESENT( ps_i ) ) THEN 
     1183            ! 
     1184            ALLOCATE( z1_ai(idim), z1_vi(idim), z1_vs(idim), ztmp(idim) ) 
     1185            ! 
     1186            WHERE( SUM( pa_i(:,:), dim=2 ) /= 0._wp )               ;   z1_ai(:) = 1._wp / SUM( pa_i(:,:), dim=2 ) 
     1187            ELSEWHERE                                               ;   z1_ai(:) = 0._wp 
     1188            END WHERE 
     1189            WHERE( SUM( pa_i(:,:) * ph_i(:,:), dim=2 ) /= 0._wp )   ;   z1_vi(:) = 1._wp / SUM( pa_i(:,:) * ph_i(:,:), dim=2 ) 
     1190            ELSEWHERE                                               ;   z1_vi(:) = 0._wp 
     1191            END WHERE 
     1192            WHERE( SUM( pa_i(:,:) * ph_s(:,:), dim=2 ) /= 0._wp )   ;   z1_vs(:) = 1._wp / SUM( pa_i(:,:) * ph_s(:,:), dim=2 ) 
     1193            ELSEWHERE                                               ;   z1_vs(:) = 0._wp 
     1194            END WHERE 
     1195            ! 
     1196            ! fill all the categories with the same value 
     1197            IF( PRESENT( pt_i  ) ) THEN 
     1198               ztmp(:) = SUM( ptmi (:,:) * pati(:,:) * phti(:,:), dim=2 ) * z1_vi(:) 
     1199               DO jl = 1, jpl 
     1200                  pt_i (:,jl) = ztmp(:) 
     1201               END DO 
     1202            ENDIF 
     1203            IF( PRESENT( pt_s  ) ) THEN 
     1204               ztmp(:) =  SUM( ptms (:,:) * pati(:,:) * phts(:,:), dim=2 ) * z1_vs(:) 
     1205               DO jl = 1, jpl 
     1206                  pt_s (:,jl) = ztmp(:) 
     1207               END DO 
     1208            ENDIF 
     1209            IF( PRESENT( pt_su ) ) THEN 
     1210               ztmp(:) = SUM( ptmsu(:,:) * pati(:,:)            , dim=2 ) * z1_ai(:) 
     1211               DO jl = 1, jpl 
     1212                  pt_su(:,jl) = ztmp(:) 
     1213               END DO 
     1214            ENDIF 
     1215            IF( PRESENT( ps_i  ) ) THEN 
     1216               ztmp(:) = SUM( psmi (:,:) * pati(:,:) * phti(:,:), dim=2 ) * z1_vi(:) 
     1217               DO jl = 1, jpl 
     1218                  ps_i (:,jl) = ztmp(:) 
     1219               END DO 
     1220            ENDIF 
     1221            ! 
     1222            DEALLOCATE( z1_ai, z1_vi, z1_vs, ztmp ) 
     1223            ! 
     1224         ENDIF 
    11011225         ! 
    11021226      ENDIF 
Note: See TracChangeset for help on using the changeset viewer.