Changeset 2788 for LMDZ5/trunk/libf/misc


Ignore:
Timestamp:
Feb 2, 2017, 7:01:50 PM (7 years ago)
Author:
dcugnet
Message:

Changes in ce0l about the way ozone forcing files are generated:

1) 3D raw input files "climoz.nc" are now handled.
2) Default behaviour is now to let the gcm interpolate in time online.

This helps to avoid huge forcing files (in particular for 3D fields).
In this case, the output files "climoz_LMDZ.nc" all have 14 records:

  • records 2-13 are obtained with records 1-12 of "climoz.nc".
  • records 1 and 14 are obtained respectively with:
    • record 12 of "climoz_m.nc" if available, of "climoz.nc" otherwise.
    • record 1 of "climoz_p.nc" if available, of "climoz.nc" otherwise.

3) If ok_daily_climoz key is TRUE, the time interpolation (one record

a day) is forced, using the 14 records described below.
This now depends on the calendar (it was on a 360 days basis only).

Changes in the gcm about the way zone forcing files are read/interpolated:

1) 3D horizontally interpolated "climoz_LMDZ.nc" files are now handled.
2) Daily files (already interpolated in time) are still handled, but their

number of records must match the expected number of days, that depends
on the calendar (records step is no longer 1/360 year).

3) 14 records monthly files are now handled (and prefered). This reduces

the I/O to a minimum and the aditional computational cost is low (simple
online linear time interpolation).

4) If adjust_tropopause key is TRUE, the input fields are stretched using

following method:

  • LMDZ dynamical tropopause is detected: Ptrop_lmdz = MAX ( P(Potential Vorticity==2PVU), P(theta==380K) )
  • file chemical tropopause is detected: Ptrop_file = P( tro3 == o3t ), where:

o3t = 91. + 28. * SIN(PI*(month-2)/6) (ppbV)

This formula comes from Thouret & al., ACP 6, 1033-1051, 2006.
The second term of the expression is multiplied by TANH(lat_deg/20.)
to account for latitude dependency.

  • File profile is streched in a +/- 5kms zone around the mean tropopause to ensure resulting tropopause matches the one of LMDZ. See procedure regr_pr_time_av for more details.
Location:
LMDZ5/trunk/libf/misc
Files:
2 added
1 edited

Legend:

Unmodified
Added
Removed

  • LMDZ5/trunk/libf/misc/slopes_m.F90

    r2440 r2788  
    1 module slopes_m
     1MODULE slopes_m
    22
    33  ! Author: Lionel GUEZ
    4 
    5   implicit none
    6 
    7   interface slopes
    8      ! This generic function computes second order slopes with Van
    9      ! Leer slope-limiting, given cell averages. Reference: Dukowicz,
    10      ! 1987, SIAM Journal on Scientific and Statistical Computing, 8,
    11      ! 305.
    12 
    13      ! The only difference between the specific functions is the rank
    14      ! of the first argument and the equal rank of the result.
    15 
    16      ! real, intent(in), rank >= 1:: f ! (n, ...) cell averages, n must be >= 1
    17      ! real, intent(in):: x(:) ! (n + 1) cell edges
    18      ! real slopes, same shape as f ! (n, ...)
    19 
    20      module procedure slopes1, slopes2, slopes3, slopes4
    21   end interface
    22 
    23   private
    24   public slopes
    25 
    26 contains
    27 
    28   pure function slopes1(f, x)
    29 
    30     real, intent(in):: f(:)
    31     real, intent(in):: x(:)
    32     real slopes1(size(f))
    33 
    34     ! Local:
    35     integer n, i
    36     real xc(size(f)) ! (n) cell centers
    37     real h
    38 
    39     !------------------------------------------------------
    40 
    41     n = size(f)
    42     forall (i = 1:n) xc(i) = (x(i) + x(i + 1)) / 2.
    43     slopes1(1) = 0.
    44     slopes1(n) = 0.
    45 
    46     do i = 2, n - 1
    47        if (f(i) >= max(f(i - 1), f(i + 1)) .or. f(i) &
    48             <= min(f(i - 1), f(i + 1))) then
    49           ! Local extremum
    50           slopes1(i) = 0.
    51        else
    52           ! (f(i - 1), f(i), f(i + 1)) strictly monotonous
    53 
    54           ! Second order slope:
    55           slopes1(i) = (f(i + 1) - f(i - 1)) / (xc(i + 1) - xc(i - 1))
    56 
    57           ! Slope limitation:
    58           h = abs(x(i + 1) - xc(i))
    59           slopes1(i) = sign(min(abs(slopes1(i)), abs(f(i + 1) - f(i)) / h, &
    60                abs(f(i) - f(i - 1)) / h), slopes1(i))
    61        end if
    62     end do
    63 
    64   end function slopes1
    65 
    66   !*************************************************************
    67 
    68   pure function slopes2(f, x)
    69 
    70     real, intent(in):: f(:, :)
    71     real, intent(in):: x(:)
    72     real slopes2(size(f, 1), size(f, 2))
    73 
    74     ! Local:
    75     integer n, i, j
    76     real xc(size(f, 1)) ! (n) cell centers
    77     real h(2:size(f, 1) - 1), delta_xc(2:size(f, 1) - 1) ! (2:n - 1)
    78 
    79     !------------------------------------------------------
    80 
    81     n = size(f, 1)
    82     forall (i = 1:n) xc(i) = (x(i) + x(i + 1)) / 2.
    83 
    84     forall (i = 2:n - 1)
    85        h(i) = abs(x(i + 1) - xc(i))
    86        delta_xc(i) = xc(i + 1) - xc(i - 1)
    87     end forall
    88 
    89     do j = 1, size(f, 2)
    90        slopes2(1, j) = 0.
    91 
    92        do i = 2, n - 1
    93           if (f(i, j) >= max(f(i - 1, j), f(i + 1, j)) .or. &
    94                f(i, j) <= min(f(i - 1, j), f(i + 1, j))) then
    95              ! Local extremum
    96              slopes2(i, j) = 0.
    97           else
    98              ! (f(i - 1, j), f(i, j), f(i + 1, j))
    99              ! strictly monotonous
    100 
    101              ! Second order slope:
    102              slopes2(i, j) = (f(i + 1, j) - f(i - 1, j)) / delta_xc(i)
    103 
    104              ! Slope limitation:
    105              slopes2(i, j) = sign(min(abs(slopes2(i, j)), &
    106                   abs(f(i + 1, j) - f(i, j)) / h(i), &
    107                   abs(f(i, j) - f(i - 1, j)) / h(i)), slopes2(i, j))
    108           end if
    109        end do
    110 
    111        slopes2(n, j) = 0.
    112     end do
    113 
    114   end function slopes2
    115 
    116   !*************************************************************
    117 
    118   pure function slopes3(f, x)
    119 
    120     real, intent(in):: f(:, :, :)
    121     real, intent(in):: x(:)
    122     real slopes3(size(f, 1), size(f, 2), size(f, 3))
    123 
    124     ! Local:
    125     integer n, i, j, k
    126     real xc(size(f, 1)) ! (n) cell centers
    127     real h(2:size(f, 1) - 1), delta_xc(2:size(f, 1) - 1) ! (2:n - 1)
    128 
    129     !------------------------------------------------------
    130 
    131     n = size(f, 1)
    132     forall (i = 1:n) xc(i) = (x(i) + x(i + 1)) / 2.
    133 
    134     forall (i = 2:n - 1)
    135        h(i) = abs(x(i + 1) - xc(i))
    136        delta_xc(i) = xc(i + 1) - xc(i - 1)
    137     end forall
    138 
    139     do k = 1, size(f, 3)
    140        do j = 1, size(f, 2)
    141           slopes3(1, j, k) = 0.
    142 
    143           do i = 2, n - 1
    144              if (f(i, j, k) >= max(f(i - 1, j, k), f(i + 1, j, k)) .or. &
    145                   f(i, j, k) <= min(f(i - 1, j, k), f(i + 1, j, k))) then
    146                 ! Local extremum
    147                 slopes3(i, j, k) = 0.
    148              else
    149                 ! (f(i - 1, j, k), f(i, j, k), f(i + 1, j, k))
    150                 ! strictly monotonous
    151 
    152                 ! Second order slope:
    153                 slopes3(i, j, k) = (f(i + 1, j, k) - f(i - 1, j, k)) &
    154                      / delta_xc(i)
    155 
    156                 ! Slope limitation:
    157                 slopes3(i, j, k) = sign(min(abs(slopes3(i, j, k)), &
    158                      abs(f(i + 1, j, k) - f(i, j, k)) / h(i), &
    159                      abs(f(i, j, k) - f(i - 1, j, k)) / h(i)), slopes3(i, j, k))
    160              end if
    161           end do
    162 
    163           slopes3(n, j, k) = 0.
    164        end do
    165     end do
    166 
    167   end function slopes3
    168 
    169   !*************************************************************
    170 
    171   pure function slopes4(f, x)
    172 
    173     real, intent(in):: f(:, :, :, :)
    174     real, intent(in):: x(:)
    175     real slopes4(size(f, 1), size(f, 2), size(f, 3), size(f, 4))
    176 
    177     ! Local:
    178     integer n, i, j, k, l
    179     real xc(size(f, 1)) ! (n) cell centers
    180     real h(2:size(f, 1) - 1), delta_xc(2:size(f, 1) - 1) ! (2:n - 1)
    181 
    182     !------------------------------------------------------
    183 
    184     n = size(f, 1)
    185     forall (i = 1:n) xc(i) = (x(i) + x(i + 1)) / 2.
    186 
    187     forall (i = 2:n - 1)
    188        h(i) = abs(x(i + 1) - xc(i))
    189        delta_xc(i) = xc(i + 1) - xc(i - 1)
    190     end forall
    191 
    192     do l = 1, size(f, 4)
    193        do k = 1, size(f, 3)
    194           do j = 1, size(f, 2)
    195              slopes4(1, j, k, l) = 0.
    196 
    197              do i = 2, n - 1
    198                 if (f(i, j, k, l) >= max(f(i - 1, j, k, l), f(i + 1, j, k, l)) &
    199                      .or. f(i, j, k, l) &
    200                      <= min(f(i - 1, j, k, l), f(i + 1, j, k, l))) then
    201                    ! Local extremum
    202                    slopes4(i, j, k, l) = 0.
    203                 else
    204                    ! (f(i - 1, j, k, l), f(i, j, k, l), f(i + 1, j, k, l))
    205                    ! strictly monotonous
    206 
    207                    ! Second order slope:
    208                    slopes4(i, j, k, l) = (f(i + 1, j, k, l) &
    209                         - f(i - 1, j, k, l)) / delta_xc(i)
    210 
    211                    ! Slope limitation:
    212                    slopes4(i, j, k, l) = sign(min(abs(slopes4(i, j, k, l)), &
    213                         abs(f(i + 1, j, k, l) - f(i, j, k, l)) / h(i), &
    214                         abs(f(i, j, k, l) - f(i - 1, j, k, l)) / h(i)), &
    215                         slopes4(i, j, k, l))
    216                 end if
    217              end do
    218 
    219              slopes4(n, j, k, l) = 0.
    220           end do
    221        end do
    222     end do
    223 
    224   end function slopes4
    225 
    226 end module slopes_m
     4  ! Extension / factorisation: David CUGNET
     5
     6  IMPLICIT NONE
     7
     8  ! Those generic function computes second order slopes with Van
     9  ! Leer slope-limiting, given cell averages. Reference: Dukowicz,
     10  ! 1987, SIAM Journal on Scientific and Statistical Computing, 8,
     11  ! 305.
     12
     13  ! The only difference between the specific functions is the rank
     14  ! of the first argument and the equal rank of the result.
     15
     16  ! slope(ix,...) acts on ix th dimension.
     17
     18  ! real, intent(in), rank >= 1:: f ! (n, ...) cell averages, n must be >= 1
     19  ! real, intent(in):: x(:) ! (n + 1) cell edges
     20  ! real slopes, same shape as f ! (n, ...)
     21  INTERFACE slopes
     22     MODULE procedure slopes1, slopes2, slopes3, slopes4, slopes5
     23  END INTERFACE
     24
     25  PRIVATE
     26  PUBLIC :: slopes
     27
     28CONTAINS
     29
     30!-------------------------------------------------------------------------------
     31!
     32PURE FUNCTION slopes1(ix, f, x)
     33!
     34!-------------------------------------------------------------------------------
     35! Arguments:
     36  INTEGER, INTENT(IN) :: ix
     37  REAL,    INTENT(IN) :: f(:)
     38  REAL,    INTENT(IN) :: x(:)
     39  REAL :: slopes1(SIZE(f,1))
     40!-------------------------------------------------------------------------------
     41! Local:
     42  INTEGER :: n, i, j, sta(2), sto(2)
     43  REAL :: xc(SIZE(f,1))                             ! (n) cell centers
     44  REAL :: h(2:SIZE(f,1)-1), delta_xc(2:SIZE(f,1)-1) ! (2:n-1)
     45  REAL :: fm, ff, fp, dx
     46!-------------------------------------------------------------------------------
     47  n=SIZE(f,ix)
     48  FORALL(i=1:n) xc(i)=(x(i)+x(i+1))/2.
     49  FORALL(i=2:n-1)
     50    h(i)=ABS(x(i+1)-xc(i)) ; delta_xc(i)=xc(i+1)-xc(i-1)
     51  END FORALL
     52  slopes1(:)=0.
     53  DO i=2,n-1
     54    ff=f(i); fm=f(i-1); fp=f(i+1)
     55    IF(ff>=MAX(fm,fp).OR.ff<=MIN(fm,fp)) THEN
     56      slopes1(i)=0.; CYCLE           !--- Local extremum
     57      !--- 2nd order slope ; (fm, ff, fp) strictly monotonous
     58      slopes1(i)=(fp-fm)/delta_xc(i)
     59      !--- Slope limitation
     60      slopes1(i) = SIGN(MIN(ABS(slopes1(i)), &
     61        ABS(fp-ff)/h(i),ABS(ff-fm)/h(i)),slopes1(i) )
     62     END IF
     63  END DO
     64
     65END FUNCTION slopes1
     66!
     67!-------------------------------------------------------------------------------
     68
     69
     70!-------------------------------------------------------------------------------
     71!
     72PURE FUNCTION slopes2(ix, f, x)
     73!
     74!-------------------------------------------------------------------------------
     75! Arguments:
     76  INTEGER, INTENT(IN) :: ix
     77  REAL,    INTENT(IN) :: f(:, :)
     78  REAL,    INTENT(IN) :: x(:)
     79  REAL :: slopes2(SIZE(f,1),SIZE(f,2))
     80!-------------------------------------------------------------------------------
     81! Local:
     82  INTEGER :: n, i, j, sta(2), sto(2)
     83  REAL, ALLOCATABLE :: xc(:)                        ! (n) cell centers
     84  REAL, ALLOCATABLE :: h(:), delta_xc(:)            ! (2:n-1)
     85  REAL :: fm, ff, fp, dx
     86!-------------------------------------------------------------------------------
     87  n=SIZE(f,ix); ALLOCATE(xc(n),h(2:n-1),delta_xc(2:n-1))
     88  FORALL(i=1:n) xc(i)=(x(i)+x(i+1))/2.
     89  FORALL(i=2:n-1)
     90    h(i)=ABS(x(i+1)-xc(i)) ; delta_xc(i)=xc(i+1)-xc(i-1)
     91  END FORALL
     92  slopes2(:,:)=0.
     93  sta=[1,1]; sta(ix)=2
     94  sto=SHAPE(f); sto(ix)=n-1
     95  DO j=sta(2),sto(2); IF(ix==2) dx=delta_xc(j)
     96    DO i=sta(1),sto(1); IF(ix==1) dx=delta_xc(i)
     97      ff=f(i,j)
     98      SELECT CASE(ix)
     99        CASE(1); fm=f(i-1,j); fp=f(i+1,j)
     100        CASE(2); fm=f(i,j-1); fp=f(i,j+1)
     101      END SELECT
     102      IF(ff>=MAX(fm,fp).OR.ff<=MIN(fm,fp)) THEN
     103        slopes2(i,j)=0.; CYCLE           !--- Local extremum
     104        !--- 2nd order slope ; (fm, ff, fp) strictly monotonous
     105        slopes2(i,j)=(fp-fm)/dx
     106        !--- Slope limitation
     107        slopes2(i,j) = SIGN(MIN(ABS(slopes2(i,j)), &
     108          ABS(fp-ff)/h(i),ABS(ff-fm)/h(i)),slopes2(i,j) )
     109       END IF
     110    END DO
     111  END DO
     112  DEALLOCATE(xc,h,delta_xc)
     113
     114END FUNCTION slopes2
     115!
     116!-------------------------------------------------------------------------------
     117
     118
     119!-------------------------------------------------------------------------------
     120!
     121PURE FUNCTION slopes3(ix, f, x)
     122!
     123!-------------------------------------------------------------------------------
     124! Arguments:
     125  INTEGER, INTENT(IN) :: ix
     126  REAL,    INTENT(IN) :: f(:, :, :)
     127  REAL,    INTENT(IN) :: x(:)
     128  REAL :: slopes3(SIZE(f,1),SIZE(f,2),SIZE(f,3))
     129!-------------------------------------------------------------------------------
     130! Local:
     131  INTEGER :: n, i, j, k, sta(3), sto(3)
     132  REAL, ALLOCATABLE :: xc(:)                        ! (n) cell centers
     133  REAL, ALLOCATABLE :: h(:), delta_xc(:)            ! (2:n-1)
     134  REAL :: fm, ff, fp, dx
     135!-------------------------------------------------------------------------------
     136  n=SIZE(f,ix); ALLOCATE(xc(n),h(2:n-1),delta_xc(2:n-1))
     137  FORALL(i=1:n) xc(i)=(x(i)+x(i+1))/2.
     138  FORALL(i=2:n-1)
     139    h(i)=ABS(x(i+1)-xc(i)) ; delta_xc(i)=xc(i+1)-xc(i-1)
     140  END FORALL
     141  slopes3(:,:,:)=0.
     142  sta=[1,1,1]; sta(ix)=2
     143  sto=SHAPE(f); sto(ix)=n-1
     144  DO k=sta(3),sto(3); IF(ix==3) dx=delta_xc(k)
     145    DO j=sta(2),sto(2); IF(ix==2) dx=delta_xc(j)
     146      DO i=sta(1),sto(1); IF(ix==1) dx=delta_xc(i)
     147        ff=f(i,j,k)
     148        SELECT CASE(ix)
     149          CASE(1); fm=f(i-1,j,k); fp=f(i+1,j,k)
     150          CASE(2); fm=f(i,j-1,k); fp=f(i,j+1,k)
     151          CASE(3); fm=f(i,j,k-1); fp=f(i,j,k+1)
     152        END SELECT
     153        IF(ff>=MAX(fm,fp).OR.ff<=MIN(fm,fp)) THEN
     154          slopes3(i,j,k)=0.; CYCLE           !--- Local extremum
     155          !--- 2nd order slope ; (fm, ff, fp) strictly monotonous
     156          slopes3(i,j,k)=(fp-fm)/dx
     157          !--- Slope limitation
     158          slopes3(i,j,k) = SIGN(MIN(ABS(slopes3(i,j,k)), &
     159            ABS(fp-ff)/h(i),ABS(ff-fm)/h(i)),slopes3(i,j,k) )
     160         END IF
     161      END DO
     162    END DO
     163  END DO
     164  DEALLOCATE(xc,h,delta_xc)
     165
     166END FUNCTION slopes3
     167!
     168!-------------------------------------------------------------------------------
     169
     170
     171!-------------------------------------------------------------------------------
     172!
     173PURE FUNCTION slopes4(ix, f, x)
     174!
     175!-------------------------------------------------------------------------------
     176! Arguments:
     177  INTEGER, INTENT(IN) :: ix
     178  REAL,    INTENT(IN) :: f(:, :, :, :)
     179  REAL,    INTENT(IN) :: x(:)
     180  REAL :: slopes4(SIZE(f,1),SIZE(f,2),SIZE(f,3),SIZE(f,4))
     181!-------------------------------------------------------------------------------
     182! Local:
     183  INTEGER :: n, i, j, k, l, m, sta(4), sto(4)
     184  REAL, ALLOCATABLE :: xc(:)                        ! (n) cell centers
     185  REAL, ALLOCATABLE :: h(:), delta_xc(:)            ! (2:n-1)
     186  REAL :: fm, ff, fp, dx
     187!-------------------------------------------------------------------------------
     188  n=SIZE(f,ix); ALLOCATE(xc(n),h(2:n-1),delta_xc(2:n-1))
     189  FORALL(i=1:n) xc(i)=(x(i)+x(i+1))/2.
     190  FORALL(i=2:n-1)
     191    h(i)=ABS(x(i+1)-xc(i)) ; delta_xc(i)=xc(i+1)-xc(i-1)
     192  END FORALL
     193  slopes4(:,:,:,:)=0.
     194  sta=[1,1,1,1]; sta(ix)=2
     195  sto=SHAPE(f); sto(ix)=n-1
     196  DO l=sta(4),sto(4); IF(ix==4) dx=delta_xc(l)
     197    DO k=sta(3),sto(3); IF(ix==3) dx=delta_xc(k)
     198      DO j=sta(2),sto(2); IF(ix==2) dx=delta_xc(j)
     199        DO i=sta(1),sto(1); IF(ix==1) dx=delta_xc(i)
     200          ff=f(i,j,k,l)
     201          SELECT CASE(ix)
     202            CASE(1); fm=f(i-1,j,k,l); fp=f(i+1,j,k,l)
     203            CASE(2); fm=f(i,j-1,k,l); fp=f(i,j+1,k,l)
     204            CASE(3); fm=f(i,j,k-1,l); fp=f(i,j,k+1,l)
     205            CASE(4); fm=f(i,j,k,l-1); fp=f(i,j,k,l+1)
     206          END SELECT
     207          IF(ff>=MAX(fm,fp).OR.ff<=MIN(fm,fp)) THEN
     208            slopes4(i,j,k,l)=0.; CYCLE           !--- Local extremum
     209            !--- 2nd order slope ; (fm, ff, fp) strictly monotonous
     210            slopes4(i,j,k,l)=(fp-fm)/dx
     211            !--- Slope limitation
     212            slopes4(i,j,k,l) = SIGN(MIN(ABS(slopes4(i,j,k,l)), &
     213              ABS(fp-ff)/h(i),ABS(ff-fm)/h(i)),slopes4(i,j,k,l) )
     214           END IF
     215        END DO
     216      END DO
     217    END DO
     218  END DO
     219  DEALLOCATE(xc,h,delta_xc)
     220
     221END FUNCTION slopes4
     222!
     223!-------------------------------------------------------------------------------
     224
     225
     226!-------------------------------------------------------------------------------
     227!
     228PURE FUNCTION slopes5(ix, f, x)
     229!
     230!-------------------------------------------------------------------------------
     231! Arguments:
     232  INTEGER, INTENT(IN) :: ix
     233  REAL,    INTENT(IN) :: f(:, :, :, :, :)
     234  REAL,    INTENT(IN) :: x(:)
     235  REAL :: slopes5(SIZE(f,1),SIZE(f,2),SIZE(f,3),SIZE(f,4),SIZE(f,5))
     236!-------------------------------------------------------------------------------
     237! Local:
     238  INTEGER :: n, i, j, k, l, m, sta(5), sto(5)
     239  REAL, ALLOCATABLE :: xc(:)                        ! (n) cell centers
     240  REAL, ALLOCATABLE :: h(:), delta_xc(:)            ! (2:n-1)
     241  REAL :: fm, ff, fp, dx
     242!-------------------------------------------------------------------------------
     243  n=SIZE(f,ix); ALLOCATE(xc(n),h(2:n-1),delta_xc(2:n-1))
     244  FORALL(i=1:n) xc(i)=(x(i)+x(i+1))/2.
     245  FORALL(i=2:n-1)
     246    h(i)=ABS(x(i+1)-xc(i)) ; delta_xc(i)=xc(i+1)-xc(i-1)
     247  END FORALL
     248  slopes5(:,:,:,:,:)=0.
     249  sta=[1,1,1,1,1]; sta(ix)=2
     250  sto=SHAPE(f);    sto(ix)=n-1
     251  DO m=sta(5),sto(5); IF(ix==5) dx=delta_xc(m)
     252    DO l=sta(4),sto(4); IF(ix==4) dx=delta_xc(l)
     253      DO k=sta(3),sto(3); IF(ix==3) dx=delta_xc(k)
     254        DO j=sta(2),sto(2); IF(ix==2) dx=delta_xc(j)
     255          DO i=sta(1),sto(1); IF(ix==1) dx=delta_xc(i)
     256            ff=f(i,j,k,l,m)
     257            SELECT CASE(ix)
     258              CASE(1); fm=f(i-1,j,k,l,m); fp=f(i+1,j,k,l,m)
     259              CASE(2); fm=f(i,j-1,k,l,m); fp=f(i,j+1,k,l,m)
     260              CASE(3); fm=f(i,j,k-1,l,m); fp=f(i,j,k+1,l,m)
     261              CASE(4); fm=f(i,j,k,l-1,m); fp=f(i,j,k,l+1,m)
     262              CASE(5); fm=f(i,j,k,l,m-1); fp=f(i,j,k,l,m+1)
     263            END SELECT
     264            IF(ff>=MAX(fm,fp).OR.ff<=MIN(fm,fp)) THEN
     265              slopes5(i,j,k,l,m)=0.; CYCLE           !--- Local extremum
     266              !--- 2nd order slope ; (fm, ff, fp) strictly monotonous
     267              slopes5(i,j,k,l,m)=(fp-fm)/dx
     268              !--- Slope limitation
     269              slopes5(i,j,k,l,m) = SIGN(MIN(ABS(slopes5(i,j,k,l,m)), &
     270                ABS(fp-ff)/h(i),ABS(ff-fm)/h(i)),slopes5(i,j,k,l,m) )
     271            END IF
     272          END DO
     273        END DO
     274      END DO
     275    END DO
     276  END DO
     277  DEALLOCATE(xc,h,delta_xc)
     278
     279END FUNCTION slopes5
     280!
     281!-------------------------------------------------------------------------------
     282
     283END MODULE slopes_m
Note: See TracChangeset for help on using the changeset viewer.