source: trunk/LMDZ.GENERIC/libf/phystd/surf_heat_transp_mod.F90 @ 1980

Last change on this file since 1980 was 1537, checked in by emillour, 9 years ago

Generic GCM:

  • Bug fix in surf_heat_transp_mod (introduced with revision 1529, 05/04/2016, modifications).
  • Initialize runoff in hydrol.F90.

EM
  • Property svn:executable set to *
File size: 14.6 KB
Line 
1!
2! $Header: /home/cvsroot/LMDZ4/libf/phylmd/ocean_slab_mod.F90,v 1.3 2008-02-04 16:24:28 fairhead Exp $
3!
4MODULE surf_heat_transp_mod
5
6IMPLICIT NONE 
7
8  ! Variables copied over from dyn3d dynamics:
9  REAL,SAVE,ALLOCATABLE :: fext(:)
10  REAL,SAVE,ALLOCATABLE :: unsairez(:)
11  REAL,SAVE,ALLOCATABLE :: unsaire(:)
12  REAL,SAVE,ALLOCATABLE :: cu(:)
13  REAL,SAVE,ALLOCATABLE :: cv(:)
14  REAL,SAVE,ALLOCATABLE :: cuvsurcv(:)
15  REAL,SAVE,ALLOCATABLE :: cvusurcu(:)
16  REAL,SAVE,ALLOCATABLE :: aire(:)
17  REAL,SAVE :: apoln
18  REAL,SAVE :: apols
19  REAL,SAVE,ALLOCATABLE :: aireu(:)
20  REAL,SAVE,ALLOCATABLE :: airev(:)
21
22  LOGICAL,SAVE :: alpha_var
23  LOGICAL,SAVE :: slab_upstream
24  REAL,SAVE,ALLOCATABLE :: zmasqu(:)
25  REAL,SAVE,ALLOCATABLE :: zmasqv(:)
26  REAL,SAVE,ALLOCATABLE :: unsfv(:)
27  REAL,SAVE,ALLOCATABLE :: unsev(:)
28  REAL,SAVE,ALLOCATABLE :: unsfu(:)
29  REAL,SAVE,ALLOCATABLE :: unseu(:)
30
31  ! Routines usable only by routines within this module:
32  PRIVATE :: gr_fi_dyn, gr_dyn_fi
33  ! Routines from dyn3d, valid on global dynamics grid only:
34  PRIVATE :: grad,diverg,gr_u_scal,gr_v_scal
35 
36CONTAINS
37 
38  SUBROUTINE ini_surf_heat_transp(ip1jm,ip1jmp1,unsairez_,fext_,unsaire_,&
39                                  cu_,cuvsurcv_,cv_,cvusurcu_, &
40                                  aire_,apoln_,apols_, &
41                                  aireu_,airev_)
42    USE mod_grid_phy_lmdz, only: nbp_lon, nbp_lat
43    IMPLICIT NONE
44    ! Transfer some variables from dyn3d dynamics
45    INTEGER,INTENT(IN) :: ip1jm
46    INTEGER,INTENT(IN) :: ip1jmp1
47    REAL,INTENT(IN) :: unsairez_(ip1jm)
48    REAL,INTENT(IN) :: fext_(ip1jm)
49    REAL,INTENT(IN) :: unsaire_(ip1jmp1)
50    REAL,INTENT(IN) :: cu_(ip1jmp1)
51    REAL,INTENT(IN) :: cuvsurcv_(ip1jm)
52    REAL,INTENT(IN) :: cv_(ip1jm)
53    REAL,INTENT(IN) :: cvusurcu_(ip1jmp1)
54    REAL,INTENT(IN) :: aire_(ip1jmp1)
55    REAL,INTENT(IN) :: apoln_
56    REAL,INTENT(IN) :: apols_
57    REAL,INTENT(IN) :: aireu_(ip1jmp1)
58    REAL,INTENT(IN) :: airev_(ip1jm)
59   
60    ! Sanity check
61    if ((ip1jm.ne.((nbp_lon+1)*(nbp_lat-1))).or. &
62        (ip1jmp1.ne.((nbp_lon+1)*nbp_lat))) then
63      write(*,*) "ini_surf_heat_transp Error: wrong array sizes"
64      stop
65    endif
66   
67    allocate(unsairez(ip1jm))
68    unsairez(:)=unsairez_(:)
69    allocate(fext(ip1jm))
70    fext(:)=fext_(:)
71    allocate(unsaire(ip1jmp1))
72    unsaire(:)=unsaire_(:)
73    allocate(cu(ip1jmp1))
74    cu(:)=cu_(:)
75    allocate(cuvsurcv(ip1jm))
76    cuvsurcv(:)=cuvsurcv_(:)
77    allocate(cv(ip1jm))
78    cv(:)=cv_(:)
79    allocate(cvusurcu(ip1jmp1))
80    cvusurcu(:)=cvusurcu_(:)
81    allocate(aire(ip1jmp1))
82    aire(:)=aire_(:)
83    apoln=apoln_
84    apols=apols_
85    allocate(aireu(ip1jmp1))
86    aireu(:)=aireu_(:)
87    allocate(airev(ip1jm))
88    airev(:)=airev_(:)
89   
90  END SUBROUTINE ini_surf_heat_transp
91 
92  SUBROUTINE ini_surf_heat_transp_mod
93    USE mod_grid_phy_lmdz, only: nbp_lon, nbp_lat
94    IMPLICIT NONE
95    INTEGER :: ip1jm, ip1jmp1
96   
97    ip1jm=(nbp_lon+1)*(nbp_lat-1)
98    ip1jmp1=(nbp_lon+1)*nbp_lat
99 
100    allocate(zmasqu(ip1jmp1))
101    allocate(zmasqv(ip1jm))
102    allocate(unsfv(ip1jm))
103    allocate(unsev(ip1jm))
104    allocate(unsfu(ip1jmp1))
105    allocate(unseu(ip1jmp1))
106
107  END SUBROUTINE ini_surf_heat_transp_mod
108
109  SUBROUTINE divgrad_phy(ngrid,nlevs,temp,delta)
110
111      USE mod_grid_phy_lmdz, ONLY: nbp_lon, nbp_lat
112      IMPLICIT NONE
113
114      INTEGER,INTENT(IN) :: ngrid, nlevs
115      REAL,INTENT(IN) :: temp(ngrid,nlevs)
116      REAL,INTENT(OUT) :: delta(ngrid,nlevs)
117      REAL delta_2d((nbp_lon+1)*nbp_lat,nlevs)
118      INTEGER :: ll
119      REAL ghx((nbp_lon+1)*nbp_lat,nlevs)
120      REAL ghy((nbp_lon+1)*(nbp_lat-1),nlevs)
121      INTEGER :: iip1,jjp1
122     
123      iip1=nbp_lon+1
124      jjp1=nbp_lat
125
126      CALL gr_fi_dyn(nlevs,ngrid,iip1,jjp1,temp,delta_2d)
127      CALL grad(nlevs,delta_2d,ghx,ghy)
128      DO ll=1,nlevs
129          ghx(:,ll)=ghx(:,ll)*zmasqu
130! pas de diffusion zonale 
131!          ghx(:,ll)=0.
132          ghy(:,ll)=ghy(:,ll)*zmasqv
133      END DO
134      CALL diverg(nlevs,ghx,ghy,delta_2d)
135      CALL gr_dyn_fi(nlevs,iip1,jjp1,ngrid,delta_2d,delta)
136
137
138  END SUBROUTINE divgrad_phy
139
140
141
142  SUBROUTINE init_masquv(ngrid,zmasq)
143
144      USE mod_grid_phy_lmdz, only: nbp_lon, nbp_lat
145      IMPLICIT NONE
146
147
148      INTEGER,INTENT(IN) :: ngrid
149      REAL zmasq(ngrid)
150      REAL zmasq_2d((nbp_lon+1)*nbp_lat)
151      REAL ff((nbp_lon+1)*(nbp_lat-1))
152      REAL eps
153      INTEGER i
154      INTEGER :: iim,iip1,jjp1,ip1jm,ip1jmp1
155     
156      iim=nbp_lon
157      iip1=nbp_lon+1
158      jjp1=nbp_lat
159      ip1jm=(nbp_lon+1)*(nbp_lat-1)
160      ip1jmp1=(nbp_lon+1)*nbp_lat
161
162! Masques u,v
163      zmasqu=1.
164      zmasqv=1.
165
166      CALL gr_fi_dyn(1,ngrid,iip1,jjp1,zmasq,zmasq_2d)
167
168      DO i=1,ip1jmp1-1
169        IF (zmasq_2d(i).GT.1e-5 .OR. zmasq_2d(i+1).GT.1e-5) THEN
170                zmasqu(i)=0.
171        ENDIF
172      END DO
173      DO i=iip1,ip1jmp1,iip1
174        zmasqu(i)=zmasqu(i-iim)
175      END DO
176      DO i=1,ip1jm
177        IF (zmasq_2d(i).GT.1e-5 .OR. zmasq_2d(i+iip1).GT.1e-5) THEN 
178                zmasqv(i)=0.
179        END IF
180      END DO
181
182
183! Coriolis (pour Ekman transp.)
184      eps=1e-5
185!      CALL getin('slab_eps',eps)
186!      print *,'epsilon=',eps
187      ff=fext*unsairez       
188      DO i=1,ip1jm
189        unsev(i)=eps/(ff(i)*ff(i)+eps**2)
190        unsfv(i)=ff(i)/(ff(i)*ff(i)+eps**2)
191      ENDDO
192      CALL gr_v_scal(1,unsfv,unsfu)
193      CALL gr_v_scal(1,unsev,unseu)
194! Alpha variable?
195!      alpha_var=.FALSE.
196!      CALL getin('slab_alphav',alpha_var)
197
198
199     
200  END SUBROUTINE init_masquv
201
202
203
204  SUBROUTINE slab_ekman2(ngrid,tx_phy,ty_phy,ts_phy,dt_phy)
205 
206      USE mod_grid_phy_lmdz, ONLY: nbp_lon, nbp_lat
207      USE slab_ice_h, ONLY: noceanmx
208
209      IMPLICIT NONE
210     
211      INTEGER,INTENT(IN) :: ngrid
212      INTEGER ij
213      REAL txv((nbp_lon+1)*(nbp_lat-1)),fluxm((nbp_lon+1)*(nbp_lat-1))
214      REAL tyv((nbp_lon+1)*(nbp_lat-1))
215      REAL fluxtm((nbp_lon+1)*(nbp_lat-1),noceanmx)
216      REAL fluxtz((nbp_lon+1)*nbp_lat,noceanmx)
217      REAL tyu((nbp_lon+1)*nbp_lat),txu((nbp_lon+1)*nbp_lat)
218      REAL fluxz((nbp_lon+1)*nbp_lat),fluxv((nbp_lon+1)*nbp_lat)
219      REAL dt((nbp_lon+1)*nbp_lat,noceanmx),ts((nbp_lon+1)*nbp_lat,noceanmx)
220      REAL tx_phy(ngrid),ty_phy(ngrid)
221      REAL dt_phy(ngrid,noceanmx),ts_phy(ngrid,noceanmx)
222
223      INTEGER iim,iip1,iip2,jjp1,ip1jm,ip1jmi1,ip1jmp1
224
225      iim=nbp_lon
226      iip1=nbp_lon+1
227      iip2=nbp_lon+2
228      jjp1=nbp_lat
229      ip1jm=(nbp_lon+1)*(nbp_lat-1) ! = iip1*jjm
230      ip1jmi1=(nbp_lon+1)*(nbp_lat-1)-(nbp_lon+1) ! = ip1jm - iip1
231      ip1jmp1=(nbp_lon+1)*nbp_lat ! = iip1*jjp1
232
233! Passage taux,y sur grilles 2d
234      CALL gr_fi_dyn(1,ngrid,iip1,jjp1,tx_phy,txu)
235      CALL gr_fi_dyn(1,ngrid,iip1,jjp1,ty_phy,tyu)
236! Passage grille u,v
237! Multiplication par f ou eps.
238      CALL gr_v_scal(1,txu,txv)
239      CALL gr_v_scal(1,tyu,tyv)
240      fluxm=tyv*unsev-txv*unsfv
241!      fluxm=-txv*unsfv
242      CALL gr_u_scal(1,txu,txu)
243      CALL gr_u_scal(1,tyu,tyu)
244      fluxz=tyu*unsfu+txu*unseu
245!      fluxz=tyu*unsfu
246           
247! Calcul de T, Tdeep
248      CALL gr_fi_dyn(2,ngrid,iip1,jjp1,ts_phy,ts)
249       
250! Flux de masse
251      fluxm=fluxm*cv*cuvsurcv*zmasqv
252      fluxz=fluxz*cu*cvusurcu*zmasqu
253! Flux de masse vertical
254      DO ij=iip2,ip1jm
255        fluxv(ij)=fluxz(ij)-fluxz(ij-1)-fluxm(ij)+fluxm(ij-iip1)
256      ENDDO
257      DO ij=iip1,ip1jmi1,iip1
258         fluxv(ij+1)=fluxv(ij+iip1)
259      END DO
260! Poles
261      fluxv(1)=-SUM(fluxm(1:iim))     
262      fluxv(ip1jmp1)=SUM(fluxm(ip1jm-iim:ip1jm-1))
263      fluxv=fluxv
264
265!Calcul flux de chaleur méridiens
266      DO ij=1,ip1jm
267          fluxtm(ij,1)=fluxm(ij)*(ts(ij+iip1,1)+ts(ij,1))/2.
268          fluxtm(ij,2)=-fluxm(ij)*(ts(ij+iip1,2)+ts(ij,2))/2.
269      END DO
270
271!Calcul flux chaleur zonaux
272      DO ij=iip2,ip1jm
273      IF (fluxz(ij).GE.0.) THEN
274             fluxtz(ij,1)=fluxz(ij)*ts(ij,1)
275             fluxtz(ij,2)=-fluxz(ij)*ts(ij+1,2)
276      ELSE
277             fluxtz(ij,1)=fluxz(ij)*ts(ij+1,1)
278             fluxtz(ij,2)=-fluxz(ij)*ts(ij,2)
279      ENDIF
280      END DO
281      DO ij=iip1*2,ip1jmp1,iip1
282             fluxtz(ij,:)=fluxtz(ij-iim,:)
283      END DO
284                   
285! Calcul de dT
286      DO ij=iip2,ip1jm
287         dt(ij,:)=fluxtz(ij-1,:)-fluxtz(ij,:)   &
288                  +fluxtm(ij,:)-fluxtm(ij-iip1,:)
289         IF (fluxv(ij).GT.0.) THEN
290           dt(ij,1)=dt(ij,1)+fluxv(ij)*ts(ij,2)
291           dt(ij,2)=dt(ij,2)-fluxv(ij)*ts(ij,2)
292         ELSE
293           dt(ij,1)=dt(ij,1)+fluxv(ij)*ts(ij,1)
294           dt(ij,2)=dt(ij,2)-fluxv(ij)*ts(ij,1)
295         ENDIF
296         dt(ij,:)=dt(ij,:)/aire(ij)
297      END DO
298      DO ij=iip1,ip1jmi1,iip1
299         dt(ij+1,:)=dt(ij+iip1,:)
300      END DO
301! Pôles
302      dt(1,:)=SUM(fluxtm(1:iim,:),dim=1)
303        IF (fluxv(1).GT.0.) THEN
304          dt(1,1)=dt(1,1)+fluxv(1)*ts(1,2)
305          dt(1,2)=dt(1,2)-fluxv(1)*ts(1,2)
306        ELSE
307          dt(1,1)=dt(1,1)+fluxv(1)*ts(1,1)
308          dt(1,2)=dt(1,2)-fluxv(1)*ts(1,1)
309        ENDIF
310      dt(1,:)=dt(1,:)/apoln
311      dt(ip1jmp1,:)=-SUM(fluxtm(ip1jm-iim:ip1jm-1,:),dim=1)
312       IF (fluxv(ip1jmp1).GT.0.) THEN
313         dt(ip1jmp1,1)=dt(ip1jmp1,1)+fluxv(ip1jmp1)*ts(ip1jmp1,2)
314         dt(ip1jmp1,2)=dt(ip1jmp1,2)-fluxv(ip1jmp1)*ts(ip1jmp1,2)
315       ELSE
316         dt(ip1jmp1,1)=dt(ip1jmp1,1)+fluxv(ip1jmp1)*ts(ip1jmp1,1)
317         dt(ip1jmp1,2)=dt(ip1jmp1,2)-fluxv(ip1jmp1)*ts(ip1jmp1,1)
318       ENDIF
319      dt(ip1jmp1,:)=dt(ip1jmp1,:)/apols
320      dt(2:iip1,1)=dt(1,1)
321      dt(2:iip1,2)=dt(1,2)
322      dt(ip1jm+1:ip1jmp1,1)=dt(ip1jmp1,1)
323      dt(ip1jm+1:ip1jmp1,2)=dt(ip1jmp1,2)
324
325! Retour grille physique
326      CALL gr_dyn_fi(2,iip1,jjp1,ngrid,dt,dt_phy)
327
328
329  END SUBROUTINE slab_ekman2
330
331!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
332
333  SUBROUTINE gr_fi_dyn(nfield,ngrid,im,jm,pfi,pdyn)
334  ! Transfer a variable on global "physics" grid to global "dynamics" grid
335  IMPLICIT NONE
336
337  INTEGER,INTENT(IN) :: im,jm,ngrid,nfield
338  REAL,INTENT(IN) :: pfi(ngrid,nfield)
339  REAL,INTENT(OUT) :: pdyn(im,jm,nfield)
340 
341  INTEGER :: i,j,ifield,ig
342 
343  DO ifield=1,nfield
344    ! Handle poles
345    DO i=1,im
346      pdyn(i,1,ifield)=pfi(1,ifield)
347      pdyn(i,jm,ifield)=pfi(ngrid,ifield)
348    ENDDO
349    ! Other points
350    DO j=2,jm-1
351      ig=2+(j-2)*(im-1)
352      CALL SCOPY(im-1,pfi(ig,ifield),1,pdyn(1,j,ifield),1)
353      pdyn(im,j,ifield)=pdyn(1,j,ifield)
354    ENDDO
355  ENDDO ! of DO ifield=1,nfield
356 
357  END SUBROUTINE gr_fi_dyn
358 
359
360
361  SUBROUTINE gr_dyn_fi(nfield,im,jm,ngrid,pdyn,pfi)
362  ! Transfer a variable on global "dynamics" grid to global "physics" grid
363  IMPLICIT NONE
364 
365  INTEGER,INTENT(IN) :: im,jm,ngrid,nfield
366  REAL,INTENT(IN) :: pdyn(im,jm,nfield)
367  REAL,INTENT(OUT) :: pfi(ngrid,nfield)
368 
369  INTEGER j,ifield,ig
370
371  ! Sanity check:
372  IF(ngrid.NE.2+(jm-2)*(im-1)) THEN
373    WRITE(*,*) "gr_dyn_fi error, wrong sizes"
374    STOP
375  ENDIF
376 
377  ! Handle poles
378  CALL SCOPY(nfield,pdyn,im*jm,pfi,ngrid)
379  CALL SCOPY(nfield,pdyn(1,jm,1),im*jm,pfi(ngrid,1),ngrid)
380  ! Other points
381  DO ifield=1,nfield
382    DO j=2,jm-1
383      ig=2+(j-2)*(im-1)
384      CALL SCOPY(im-1,pdyn(1,j,ifield),1,pfi(ig,ifield),1)
385    ENDDO
386  ENDDO
387
388  END SUBROUTINE gr_dyn_fi
389
390
391
392  SUBROUTINE  grad(klevel,pg,pgx,pgy)
393  ! compute the covariant components x,y of the gradient of pg
394  USE mod_grid_phy_lmdz, ONLY: nbp_lon, nbp_lat
395  IMPLICIT NONE
396 
397  INTEGER,INTENT(IN) :: klevel
398  REAL,INTENT(IN) :: pg((nbp_lon+1)*nbp_lat,klevel)
399  REAL,INTENT(OUT) :: pgx((nbp_lon+1)*nbp_lat,klevel)
400  REAL,INTENT(OUT) :: pgy((nbp_lon+1)*(nbp_lat-1),klevel)
401
402  INTEGER :: l,ij
403  INTEGER :: iim,iip1,ip1jm,ip1jmp1
404 
405  iim=nbp_lon
406  iip1=nbp_lon+1
407  ip1jm=(nbp_lon+1)*(nbp_lat-1) ! = iip1*jjm
408  ip1jmp1=(nbp_lon+1)*nbp_lat ! = iip1*jjp1
409 
410  DO l=1,klevel
411    DO ij=1,ip1jmp1-1
412      pgx(ij,l)=pg(ij+1,l)-pg(ij,l)
413    ENDDO
414    ! correction for pgx(ip1,j,l) ...
415    ! ... pgx(iip1,j,l)=pgx(1,j,l) ...
416    DO ij=iip1,ip1jmp1,iip1
417      pgx(ij,l)=pgx(ij-iim,l)
418    ENDDO
419    DO ij=1,ip1jm
420      pgy(ij,l)=pg(ij,l)-pg(ij+iip1,l)
421    ENDDO
422  ENDDO
423 
424  END SUBROUTINE grad
425
426
427
428  SUBROUTINE diverg(klevel,x,y,div)
429  ! compute the divergence of a vector field of components
430  ! x,y. y and y being covriant components
431  USE mod_grid_phy_lmdz, ONLY: nbp_lon, nbp_lat
432  IMPLICIT NONE
433 
434  INTEGER,INTENT(IN) :: klevel
435  REAL,INTENT(IN) :: x((nbp_lon+1)*nbp_lat,klevel)
436  REAL,INTENT(IN) :: y((nbp_lon+1)*(nbp_lat-1),klevel)
437  REAL,INTENT(OUT) :: div((nbp_lon+1)*nbp_lat,klevel)
438 
439  INTEGER :: l,ij
440  INTEGER :: iim,iip1,iip2,ip1jm,ip1jmp1,ip1jmi1
441 
442  REAL :: aiy1(nbp_lon+1),aiy2(nbp_lon+1)
443  REAL :: sumypn,sumyps
444  REAL,EXTERNAL :: SSUM
445 
446  iim=nbp_lon
447  iip1=nbp_lon+1
448  iip2=nbp_lon+2
449  ip1jm=(nbp_lon+1)*(nbp_lat-1) ! = iip1*jjm
450  ip1jmp1=(nbp_lon+1)*nbp_lat ! = iip1*jjp1
451  ip1jmi1=(nbp_lon+1)*(nbp_lat-1)-(nbp_lon+1) ! = ip1jm - iip1
452 
453  DO l=1,klevel
454    DO ij=iip2,ip1jm-1
455      div(ij+1,l)= &
456        cvusurcu(ij+1)*x(ij+1,l)-cvusurcu(ij)*x(ij,l)+ &
457        cuvsurcv(ij-iim)*y(ij-iim,l)-cuvsurcv(ij+1)*y(ij+1,l)
458    ENDDO
459    ! correction for div(1,j,l) ...
460    ! ... div(1,j,l)= div(iip1,j,l) ...
461    DO ij=iip2,ip1jm,iip1
462      div(ij,l)=div(ij+iim,l)
463    ENDDO
464    ! at the poles
465    DO ij=1,iim
466      aiy1(ij)=cuvsurcv(ij)*y(ij,l)
467      aiy2(ij)=cuvsurcv(ij+ip1jmi1)*y(ij+ip1jmi1,l)
468    ENDDO
469    sumypn=SSUM(iim,aiy1,1)/apoln
470    sumyps=SSUM(iim,aiy2,1)/apols
471    DO ij=1,iip1
472      div(ij,l)=-sumypn
473      div(ij+ip1jm,l)=sumyps
474    ENDDO
475  ENDDO ! of DO l=1,klevel
476 
477  !!! CALL filtreg( div, jjp1, klevel, 2, 2, .TRUE., 1 )
478 
479  DO l=1,klevel
480    DO ij=iip2,ip1jm
481      div(ij,l)=div(ij,l)*unsaire(ij)
482    ENDDO
483  ENDDO
484 
485  END SUBROUTINE diverg
486
487
488
489  SUBROUTINE gr_u_scal(nx,x_u,x_scal)
490  USE mod_grid_phy_lmdz, ONLY: nbp_lon, nbp_lat
491  IMPLICIT NONE
492 
493  INTEGER,INTENT(IN) :: nx
494  REAL,INTENT(IN) :: x_u((nbp_lon+1)*nbp_lat,nx)
495  REAL,INTENT(OUT) :: x_scal((nbp_lon+1)*nbp_lat,nx)
496 
497  INTEGER :: l,ij
498  INTEGER :: iip1,jjp1,ip1jmp1
499 
500  iip1=nbp_lon+1
501  jjp1=nbp_lat
502  ip1jmp1=(nbp_lon+1)*nbp_lat ! = iip1*jjp1
503 
504  DO l=1,nx
505    DO ij=ip1jmp1,2,-1
506      x_scal(ij,l)= &
507                   (aireu(ij)*x_u(ij,l)+aireu(ij-1)*x_u(ij-1,l)) &
508                  /(aireu(ij)+aireu(ij-1))
509    ENDDO
510  ENDDO
511 
512  CALL SCOPY(nx*jjp1,x_scal(iip1,1),iip1,x_scal(1,1),iip1)
513
514  END SUBROUTINE gr_u_scal
515
516
517
518  SUBROUTINE gr_v_scal(nx,x_v,x_scal)
519  USE mod_grid_phy_lmdz, ONLY: nbp_lon, nbp_lat
520  IMPLICIT NONE
521 
522  INTEGER,INTENT(IN) :: nx
523  REAL,INTENT(IN) :: x_v((nbp_lon+1)*(nbp_lat-1),nx)
524  REAL,INTENT(OUT) :: x_scal((nbp_lon+1)*nbp_lat,nx)
525 
526  INTEGER :: l,ij
527  INTEGER :: iip1,iip2,ip1jm,ip1jmp1
528 
529  iip1=nbp_lon+1
530  iip2=nbp_lon+2
531  ip1jm=(nbp_lon+1)*(nbp_lat-1) ! = iip1*jjm
532  ip1jmp1=(nbp_lon+1)*nbp_lat ! = iip1*jjp1
533 
534  DO l=1,nx
535    DO ij=iip2,ip1jm
536      x_scal(ij,l)= &
537                   (airev(ij-iip1)*x_v(ij-iip1,l)+airev(ij)*x_v(ij,l)) &
538                  /(airev(ij-iip1)+airev(ij))
539    ENDDO
540    DO ij=1,iip1
541      x_scal(ij,l)=0.
542    ENDDO
543    DO ij=ip1jm+1,ip1jmp1
544      x_scal(ij,l)=0.
545    ENDDO
546  ENDDO
547
548  END SUBROUTINE gr_v_scal
549 
550END MODULE surf_heat_transp_mod
551
552
553
Note: See TracBrowser for help on using the repository browser.