source: LMDZ4/trunk/libf/dyn3d/extrapol.F @ 5051

Last change on this file since 5051 was 1403, checked in by Laurent Fairhead, 14 years ago

Merged LMDZ4V5.0-dev branch changes r1292:r1399 to trunk.

Validation:
Validation consisted in compiling the HEAD revision of the trunk,
LMDZ4V5.0-dev branch and the merged sources and running different
configurations on local and SX8 machines comparing results.

Local machine: bench configuration, 32x24x11, gfortran

  • IPSLCM5A configuration (comparison between trunk and merged sources):
    • numerical convergence on dynamical fields over 3 days
    • start files are equivalent (except for RN and PB fields)
    • daily history files equivalent
  • MH07 configuration, new physics package (comparison between LMDZ4V5.0-dev branch and merged sources):
    • numerical convergence on dynamical fields over 3 days
    • start files are equivalent (except for RN and PB fields)
    • daily history files equivalent

SX8 machine (brodie), 96x95x39 on 4 processors:

  • IPSLCM5A configuration:
    • start files are equivalent (except for RN and PB fields)
    • monthly history files equivalent
  • MH07 configuration:
    • start files are equivalent (except for RN and PB fields)
    • monthly history files equivalent

Changes to the makegcm and create_make_gcm scripts to take into account
main programs in F90 files

Fusion de la branche LMDZ4V5.0-dev (r1292:r1399) au tronc principal

Validation:
La validation a consisté à compiler la HEAD de le trunk et de la banche
LMDZ4V5.0-dev et les sources fusionnées et de faire tourner le modéle selon
différentes configurations en local et sur SX8 et de comparer les résultats

En local: 32x24x11, config bench/gfortran

  • pour une config IPSLCM5A (comparaison tronc/fusion):
    • convergence numérique sur les champs dynamiques après 3 jours
    • restart et restartphy égaux (à part sur RN et Pb)
    • fichiers histoire égaux
  • pour une config nlle physique (MH07) (comparaison LMDZ4v5.0-dev/fusion):
    • convergence numérique sur les champs dynamiques après 3 jours
    • restart et restartphy égaux
    • fichiers histoire équivalents

Sur brodie, 96x95x39 sur 4 proc:

  • pour une config IPSLCM5A:
    • restart et restartphy égaux (à part sur RN et PB)
    • pas de différence dans les fichiers histmth.nc
  • pour une config MH07
    • restart et restartphy égaux (à part sur RN et PB)
    • pas de différence dans les fichiers histmth.nc

Changement sur makegcm et create_make-gcm pour pouvoir prendre en compte des
programmes principaux en *F90
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 5.0 KB
RevLine 
[524]1!
[1403]2! $Id: extrapol.F 1403 2010-07-01 09:02:53Z evignon $
[524]3!
4C
5C
6      SUBROUTINE extrapol (pfild, kxlon, kylat, pmask,
7     .                   norsud, ldper, knbor, pwork)
8      IMPLICIT none
9c
10c OASIS routine (Adaptation: Laurent Li, le 14 mars 1997)
11c Fill up missed values by using the neighbor points
12c
13      INTEGER kxlon, kylat ! longitude and latitude dimensions (Input)
14      INTEGER knbor ! minimum neighbor number (Input)
15      LOGICAL norsud ! True if field is from North to South (Input)
16      LOGICAL ldper ! True if take into account the periodicity (Input)
17      REAL pmask ! mask value (Input)
18      REAL pfild(kxlon,kylat) ! field to be extrapolated (Input/Output)
19      REAL pwork(kxlon,kylat) ! working space
20c
21      REAL zwmsk
22      INTEGER incre, idoit, i, j, k, inbor, ideb, ifin, ilon, jlat
23      INTEGER ix(9), jy(9) ! index arrays for the neighbors coordinates
24      REAL zmask(9)
25C
26C  We search over the eight closest neighbors
27C
28C            j+1  7  8  9
29C              j  4  5  6    Current point 5 --> (i,j)
30C            j-1  1  2  3
31C                i-1 i i+1
32c
33c
34      IF (norsud) THEN
35         DO j = 1, kylat
36         DO i = 1, kxlon
37            pwork(i,j) = pfild(i,kylat-j+1)
38         ENDDO
39         ENDDO
40         DO j = 1, kylat
41         DO i = 1, kxlon
42            pfild(i,j) = pwork(i,j)
43         ENDDO
44         ENDDO
45      ENDIF
46c
47      incre = 0
48c
49      DO j = 1, kylat
50      DO i = 1, kxlon
51         pwork(i,j) = pfild(i,j)
52      ENDDO
53      ENDDO
54c
55C* To avoid problems in floating point tests
56      zwmsk = pmask - 1.0
57c
58200   CONTINUE
59      incre = incre + 1
60      DO 99999 j = 1, kylat
61      DO 99999 i = 1, kxlon
62      IF (pfild(i,j).GT. zwmsk) THEN
63         pwork(i,j) = pfild(i,j)
64         inbor = 0
65         ideb = 1
66         ifin = 9
67C
68C* Fill up ix array
69         ix(1) = MAX (1,i-1)
70         ix(2) = i
71         ix(3) = MIN (kxlon,i+1)
72         ix(4) = MAX (1,i-1)
73         ix(5) = i
74         ix(6) = MIN (kxlon,i+1)
75         ix(7) = MAX (1,i-1)
76         ix(8) = i
77         ix(9) = MIN (kxlon,i+1)
78C
79C* Fill up iy array
80         jy(1) = MAX (1,j-1)
81         jy(2) = MAX (1,j-1)
82         jy(3) = MAX (1,j-1)
83         jy(4) = j
84         jy(5) = j
85         jy(6) = j
86         jy(7) = MIN (kylat,j+1)
87         jy(8) = MIN (kylat,j+1)
88         jy(9) = MIN (kylat,j+1)
89C
90C* Correct latitude bounds if southernmost or northernmost points
91         IF (j .EQ. 1) ideb = 4
92         IF (j .EQ. kylat) ifin = 6
93C
94C* Account for periodicity in longitude
95C
96         IF (ldper) THEN
97            IF (i .EQ. kxlon) THEN
98               ix(3) = 1
99               ix(6) = 1
100               ix(9) = 1
101            ELSE IF (i .EQ. 1) THEN
102               ix(1) = kxlon
103               ix(4) = kxlon
104               ix(7) = kxlon
105            ENDIF
106         ELSE
107            IF (i .EQ. 1) THEN
108               ix(1) = i
109               ix(2) = i + 1
110               ix(3) = i
111               ix(4) = i + 1
112               ix(5) = i
113               ix(6) = i + 1
114            ENDIF
115            IF (i .EQ. kxlon) THEN
116               ix(1) = i -1
117               ix(2) = i
118               ix(3) = i - 1
119               ix(4) = i
120               ix(5) = i - 1
121               ix(6) = i
122            ENDIF
123C
124            IF (i .EQ. 1 .OR. i .EQ. kxlon) THEN
125               jy(1) = MAX (1,j-1)
126               jy(2) = MAX (1,j-1)
127               jy(3) = j
128               jy(4) = j
129               jy(5) = MIN (kylat,j+1)
130               jy(6) = MIN (kylat,j+1)
131C
132               ideb = 1
133               ifin = 6
134               IF (j .EQ. 1) ideb = 3
135               IF (j .EQ. kylat) ifin = 4
136            ENDIF
137         ENDIF ! end for ldper test
138C
139C* Find unmasked neighbors
140C
141         DO 230 k = ideb, ifin
142            zmask(k) = 0.
143            ilon = ix(k)
144            jlat = jy(k)
145            IF (pfild(ilon,jlat) .LT. zwmsk) THEN
146               zmask(k) = 1.
147               inbor = inbor + 1
148            ENDIF
149 230     CONTINUE
150C
151C* Not enough points around point P are unmasked; interpolation on P
152C  will be done in a future call to extrap.
153C
154         IF (inbor .GE. knbor) THEN
155            pwork(i,j) = 0.
156            DO k = ideb, ifin
157               ilon = ix(k)
158               jlat = jy(k)
159               pwork(i,j) = pwork(i,j)
[1403]160     $                      + pfild(ilon,jlat) * zmask(k)/REAL(inbor)
[524]161            ENDDO
162         ENDIF
163C
164      ENDIF
16599999 CONTINUE
166C
167C*    3. Writing back unmasked field in pfild
168C        ------------------------------------
169C
170C* pfild then contains:
171C     - Values which were not masked
172C     - Interpolated values from the inbor neighbors
173C     - Values which are not yet interpolated
174C
175      idoit = 0
176      DO j = 1, kylat
177      DO i = 1, kxlon
178         IF (pwork(i,j) .GT. zwmsk) idoit = idoit + 1
179         pfild(i,j) = pwork(i,j)
180      ENDDO
181      ENDDO
182c
183      IF (idoit .ne. 0) GOTO 200
184ccc      PRINT*, "Number of extrapolation steps incre =", incre
185c
186      IF (norsud) THEN
187         DO j = 1, kylat
188         DO i = 1, kxlon
189            pwork(i,j) = pfild(i,kylat-j+1)
190         ENDDO
191         ENDDO
192         DO j = 1, kylat
193         DO i = 1, kxlon
194            pfild(i,j) = pwork(i,j)
195         ENDDO
196         ENDDO
197      ENDIF
198c
199      RETURN
200      END
Note: See TracBrowser for help on using the repository browser.