New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 13141 for NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_oce_interp.F90 – NEMO

Ignore:
Timestamp:
2020-06-22T18:27:34+02:00 (4 years ago)
Author:
jchanut
Message:

#2129, corrections/add ons to initial state interpolation with AGRIF
1) add namelist flag for child grid initial state interpolation - ice not considered yet
2) provide depths and not thicknesses as inputs to vertical linear interpolation
3) extend initial state interpolation to a restart scenario for parent grid (warning should be added in that case in order to prevent users doing this at each model restart...)
The online interpolation seems to work fine in the VORTEX case (provided 0. is not considered as a special value in the initial velocity field, i.e. ln_spc_dyn=F)

File:
1 edited

Legend:

Unmodified
Added
Removed

  • NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_oce_interp.F90

    r13026 r13141  
    717717      REAL(wp) :: zhtot 
    718718      REAL(wp), DIMENSION(k1:k2,1:jpts) :: tabin 
    719       REAL(wp), DIMENSION(k1:k2) :: h_in 
    720       REAL(wp), DIMENSION(1:jpk) :: h_out 
     719      REAL(wp), DIMENSION(k1:k2) :: h_in, z_in 
     720      REAL(wp), DIMENSION(1:jpk) :: h_out, z_out 
    721721      !!---------------------------------------------------------------------- 
    722722 
    … …  
    790790                     tabin(jk,:) = ptab(ji,jj,jk,n1:n2-1) 
    791791                  END DO 
     792                  z_in(1) = 0.5_wp * h_in(1) - zhtot + ht0_parent(ji,jj) 
     793                  DO jk=2,N_in 
     794                     z_in(jk) = z_in(jk-1) + 0.5_wp * h_in(jk) 
     795                  ENDDO 
     796 
    792797                  N_out = 0 
    793798                  DO jk=1,jpk ! jpk of child grid 
    … …  
    796801                     h_out(jk) = e3t(ji,jj,jk,Krhs_a) 
    797802                  ENDDO 
     803 
     804                  z_out(1) = 0.5_wp * h_out(1) - SUM(h_out(1:N_out)) + ht_0(ji,jj) 
     805                  DO jk=2,N_out 
     806                     z_out(jk) = z_out(jk-1) + 0.5_wp * h_out(jk) 
     807                  ENDDO 
     808 
    798809                  IF (N_in*N_out > 0) THEN 
    799810                     IF( l_ini_child ) THEN 
    800                         CALL remap_linear(tabin(1:N_in,1:jpts),h_in(1:N_in),ts(ji,jj,1:N_out,1:jpts,Krhs_a),          & 
    801                                       &   h_out(1:N_out),N_in,N_out,jpts)   
     811                        CALL remap_linear(tabin(1:N_in,1:jpts),z_in(1:N_in),ts(ji,jj,1:N_out,1:jpts,Krhs_a),          & 
     812                                      &   z_out(1:N_out),N_in,N_out,jpts)   
    802813                     ELSE  
    803814                        CALL reconstructandremap(tabin(1:N_in,1:jpts),h_in(1:N_in),ts(ji,jj,1:N_out,1:jpts,Krhs_a),   & 
    … …  
    850861      REAL(wp) :: zrhoy, zhtot 
    851862      ! vertical interpolation: 
    852       REAL(wp), DIMENSION(k1:k2) :: tabin, h_in 
    853       REAL(wp), DIMENSION(1:jpk) :: h_out 
     863      REAL(wp), DIMENSION(k1:k2) :: tabin, h_in, z_in 
     864      REAL(wp), DIMENSION(1:jpk) :: h_out, z_out 
    854865      INTEGER  :: N_in, N_out,item 
    855866      REAL(wp) :: h_diff 
    … …  
    924935                     ENDIF 
    925936                 ENDDO 
     937                 z_in(1) = 0.5_wp * h_in(1) - zhtot + hu0_parent(ji,jj)  
     938                 DO jk=2,N_in 
     939                    z_in(jk) = z_in(jk-1) + 0.5_wp * h_in(jk) 
     940                 ENDDO 
    926941                      
    927942                 N_out = 0 
    … …  
    931946                    h_out(N_out) = e3u(ji,jj,jk,Krhs_a) 
    932947                 ENDDO 
     948 
     949                 z_out(1) = 0.5_wp * h_out(1) - SUM(h_out(1:N_out)) + hu_0(ji,jj) 
     950                 DO jk=2,N_out 
     951                    z_out(jk) = z_out(jk-1) + 0.5_wp * h_out(jk)  
     952                 ENDDO   
     953 
    933954                 IF (N_in*N_out > 0) THEN 
    934955                     IF( l_ini_child ) THEN 
    935                         CALL remap_linear       (tabin(1:N_in),h_in(1:N_in),uu(ji,jj,1:N_out,Krhs_a),h_out(1:N_out),N_in,N_out,1) 
     956                        CALL remap_linear       (tabin(1:N_in),z_in(1:N_in),uu(ji,jj,1:N_out,Krhs_a),z_out(1:N_out),N_in,N_out,1) 
    936957                     ELSE 
    937958                        CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),uu(ji,jj,1:N_out,Krhs_a),h_out(1:N_out),N_in,N_out,1) 
    … …  
    964985      REAL(wp) :: zrhox 
    965986      ! vertical interpolation: 
    966       REAL(wp), DIMENSION(k1:k2) :: tabin, h_in 
    967       REAL(wp), DIMENSION(1:jpk) :: h_out 
     987      REAL(wp), DIMENSION(k1:k2) :: tabin, h_in, z_in 
     988      REAL(wp), DIMENSION(1:jpk) :: h_out, z_out 
    968989      INTEGER  :: N_in, N_out, item 
    969990      REAL(wp) :: h_diff, zhtot 
    … …  
    10291050                     ENDIF 
    10301051                     zhtot = zhtot + h_in(jk) 
    1031                     IF( h_in(jk) .GT. 0. ) THEN 
    1032                      tabin(jk) = ptab(ji,jj,jk,1)/(e1v(ji,jj)*zrhox*h_in(jk)) 
    1033                     ELSE 
    1034                      tabin(jk)  = 0. 
    1035                     ENDIF  
    1036                  ENDDO 
    1037              
     1052                     IF( h_in(jk) .GT. 0. ) THEN 
     1053                       tabin(jk) = ptab(ji,jj,jk,1)/(e1v(ji,jj)*zrhox*h_in(jk)) 
     1054                     ELSE 
     1055                       tabin(jk)  = 0. 
     1056                     ENDIF  
     1057                  ENDDO 
     1058 
     1059                  z_in(1) = 0.5_wp * h_in(1) - zhtot + hv0_parent(ji,jj) 
     1060                  DO jk=2,N_in 
     1061                     z_in(jk) = z_in(jk-1) + 0.5_wp * h_in(jk) 
     1062                  ENDDO 
     1063 
    10381064                  N_out = 0 
    10391065                  DO jk=1,jpk 
    1040                      if (vmask(ji,jj,jk) == 0) EXIT 
     1066                     IF (vmask(ji,jj,jk) == 0) EXIT 
    10411067                     N_out = N_out + 1 
    10421068                     h_out(N_out) = e3v(ji,jj,jk,Krhs_a) 
    1043                   END DO 
     1069                  ENDDO 
     1070 
     1071                  z_out(1) = 0.5_wp * h_out(1) - SUM(h_out(1:N_out)) + hv_0(ji,jj) 
     1072                  DO jk=2,N_out 
     1073                     z_out(jk) = z_out(jk-1) + 0.5_wp * h_out(jk) 
     1074                  ENDDO 
     1075  
    10441076                  IF (N_in*N_out > 0) THEN 
    10451077                     IF( l_ini_child ) THEN 
    1046                         CALL remap_linear       (tabin(1:N_in),h_in(1:N_in),vv(ji,jj,1:N_out,Krhs_a),h_out(1:N_out),N_in,N_out,1) 
     1078                        CALL remap_linear       (tabin(1:N_in),z_in(1:N_in),vv(ji,jj,1:N_out,Krhs_a),z_out(1:N_out),N_in,N_out,1) 
    10471079                     ELSE 
    10481080                        CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),vv(ji,jj,1:N_out,Krhs_a),h_out(1:N_out),N_in,N_out,1) 
Note: See TracChangeset for help on using the changeset viewer.