source: LMDZ4/trunk/libf/bibio/regr1_lint_m.F90 @ 5408

Last change on this file since 5408 was 1279, checked in by Laurent Fairhead, 15 years ago

Merged LMDZ4-dev branch changes r1241:1278 into the trunk
Running trunk and LMDZ4-dev in LMDZOR configuration on local
machine (sequential) and SX8 (4-proc) yields identical results
(restart and restartphy are identical binarily)
Log history from r1241 to r1278 is available by switching to
source:LMDZ4/branches/LMDZ4-dev-20091210

File size: 2.6 KB
Line 
1! $Id$
2module regr1_lint_m
3
4  ! Author: Lionel GUEZ
5
6  implicit none
7
8  interface regr1_lint
9     ! Each procedure regrids by linear interpolation.
10     ! The regridding operation is done on the first dimension of the
11     ! input array.
12     ! The difference betwwen the procedures is the rank of the first argument.
13     module procedure regr11_lint, regr12_lint
14  end interface
15
16  private
17  public regr1_lint
18
19contains
20
21  function regr11_lint(vs, xs, xt) result(vt)
22
23    ! "vs" has rank 1.
24
25    use assert_eq_m, only: assert_eq
26    use interpolation, only: hunt !!, polint
27
28    real, intent(in):: vs(:)
29    ! (values of the function at source points "xs")
30
31    real, intent(in):: xs(:)
32    ! (abscissas of points in source grid, in strictly monotonic order)
33
34    real, intent(in):: xt(:)
35    ! (abscissas of points in target grid)
36
37    real vt(size(xt)) ! values of the function on the target grid
38
39    ! Variables local to the procedure:
40    integer is, it, ns
41    integer is_b ! "is" bound between 1 and "ns - 1"
42
43    !--------------------------------------
44
45    ns = assert_eq(size(vs), size(xs), "regr11_lint ns")
46
47    is = -1 ! go immediately to bisection on first call to "hunt"
48
49    do it = 1, size(xt)
50       call hunt(xs, xt(it), is)
51       is_b = min(max(is, 1), ns - 1)
52!!       call polint(xs(is_b:is_b+1), vs(is_b:is_b+1), xt(it), vt(it))
53       vt(it) = ((xs(is_b+1) - xt(it)) * vs(is_b) &
54            + (xt(it) - xs(is_b)) * vs(is_b+1)) / (xs(is_b+1) - xs(is_b))
55    end do
56
57  end function regr11_lint
58
59  !*********************************************************
60
61  function regr12_lint(vs, xs, xt) result(vt)
62
63    ! "vs" has rank 2.
64
65    use assert_eq_m, only: assert_eq
66    use interpolation, only: hunt
67
68    real, intent(in):: vs(:, :)
69    ! (values of the function at source points "xs")
70
71    real, intent(in):: xs(:)
72    ! (abscissas of points in source grid, in strictly monotonic order)
73
74    real, intent(in):: xt(:)
75    ! (abscissas of points in target grid)
76
77    real vt(size(xt), size(vs, 2)) ! values of the function on the target grid
78
79    ! Variables local to the procedure:
80    integer is, it, ns
81    integer is_b ! "is" bound between 1 and "ns - 1"
82
83    !--------------------------------------
84
85    ns = assert_eq(size(vs, 1), size(xs), "regr12_lint ns")
86
87    is = -1 ! go immediately to bisection on first call to "hunt"
88
89    do it = 1, size(xt)
90       call hunt(xs, xt(it), is)
91       is_b = min(max(is, 1), ns - 1)
92       vt(it, :) = ((xs(is_b+1) - xt(it)) * vs(is_b, :) &
93            + (xt(it) - xs(is_b)) * vs(is_b+1, :)) / (xs(is_b+1) - xs(is_b))
94    end do
95
96  end function regr12_lint
97
98end module regr1_lint_m
Note: See TracBrowser for help on using the repository browser.