source: LMDZ4/branches/LMDZ4V5.0-dev/libf/dyn3d/sw_case_williamson91_6.F @ 1363

Last change on this file since 1363 was 1363, checked in by Ehouarn Millour, 14 years ago

Added possibility to run in "Shallow Water" mode. SW mode is automatically set if the number of atmospheric layers is 1.
To use SW mode, the model should be compiled without physics (makelmdz_fcm -p nophys) and/or without calls to the physics (i.e. set iflag_phys=0 in gcm.def).

-Updated some definitions and comments in run.def & gcm.def

-Fixed misplaced #ifdef CPP_EARTH in calfis.F + some write(lunout,*)

-Specific initialisation of fields for SW are done in sw_case_williamson91_6 (called by iniacademic, when read_start=.false.)

  • Checked (on Vargas & Brodie) that these changes don't alter usual bench GCM outputs.

EM
File size: 4.3 KB
Line 
1!
2! $Id $
3!
4      SUBROUTINE sw_case_williamson91_6(vcov,ucov,teta,masse,ps)
5
6c=======================================================================
7c
8c   Author:    Thomas Dubos      original: 26/01/2010
9c   -------
10c
11c   Subject:
12c   ------
13c   Realise le cas-test 6 de Williamson et al. (1991) : onde de Rossby-Haurwitz
14c
15c   Method:
16c   --------
17c
18c   Interface:
19c   ----------
20c
21c      Input:
22c      ------
23c
24c      Output:
25c      -------
26c
27c=======================================================================
28      IMPLICIT NONE
29c-----------------------------------------------------------------------
30c   Declararations:
31c   ---------------
32
33#include "dimensions.h"
34#include "paramet.h"
35#include "comvert.h"
36#include "comconst.h"
37#include "comgeom.h"
38#include "iniprint.h"
39
40c   Arguments:
41c   ----------
42
43c   variables dynamiques
44      REAL vcov(ip1jm,llm),ucov(ip1jmp1,llm) ! vents covariants
45      REAL teta(ip1jmp1,llm)                 ! temperature potentielle
46      REAL ps(ip1jmp1)                       ! pression  au sol
47      REAL masse(ip1jmp1,llm)                ! masse d'air
48      REAL phis(ip1jmp1)                     ! geopotentiel au sol
49
50c   Local:
51c   ------
52
53      REAL p (ip1jmp1,llmp1  )               ! pression aux interfac.des couches
54      REAL pks(ip1jmp1)                      ! exner au  sol
55      REAL pk(ip1jmp1,llm)                   ! exner au milieu des couches
56      REAL pkf(ip1jmp1,llm)                  ! exner filt.au milieu des couches
57      REAL alpha(ip1jmp1,llm),beta(ip1jmp1,llm)
58
59      REAL :: sinth,costh,costh2, Ath,Bth,Cth, lon,dps
60      INTEGER i,j,ij
61
62      REAL, PARAMETER    :: rho=1 ! masse volumique de l'air (arbitraire)
63      REAL, PARAMETER    :: K    = 7.848e-6  ! K = \omega
64      REAL, PARAMETER    :: gh0  = 9.80616 * 8e3
65      INTEGER, PARAMETER :: R0=4, R1=R0+1, R2=R0+2         ! mode 4
66c NB : rad = 6371220 dans W91 (6371229 dans LMDZ)
67c      omeg = 7.292e-5 dans W91 (7.2722e-5 dans LMDZ)
68 
69      IF(0==0) THEN
70c Williamson et al. (1991) : onde de Rossby-Haurwitz
71         teta = preff/rho/cpp
72c geopotentiel (pression de surface)
73         do j=1,jjp1
74            costh2 = cos(rlatu(j))**2
75            Ath = (R0+1)*(costh2**2) + (2*R0*R0-R0-2)*costh2 - 2*R0*R0
76            Ath = .25*(K**2)*(costh2**(R0-1))*Ath
77            Ath = .5*K*(2*omeg+K)*costh2 + Ath
78            Bth = (R1*R1+1)-R1*R1*costh2
79            Bth = 2*(omeg+K)*K/(R1*R2) * (costh2**(R0/2))*Bth
80            Cth = R1*costh2 - R2
81            Cth = .25*K*K*(costh2**R0)*Cth
82            do i=1,iip1
83               ij=(j-1)*iip1+i
84               lon = rlonv(i)
85               dps = Ath + Bth*cos(R0*lon) + Cth*cos(2*R0*lon)
86               ps(ij) = rho*(gh0 + (rad**2)*dps)
87            enddo
88         enddo
89         write(lunout,*) 'W91 ps', MAXVAL(ps), MINVAL(ps)
90c vitesse zonale ucov
91         do j=1,jjp1
92            costh  = cos(rlatu(j))
93            costh2 = costh**2
94            Ath = rad*K*costh
95            Bth = R0*(1-costh2)-costh2
96            Bth = rad*K*Bth*(costh**(R0-1))
97            do i=1,iip1
98               ij=(j-1)*iip1+i
99               lon = rlonu(i)
100               ucov(ij,1) = (Ath + Bth*cos(R0*lon))
101            enddo
102         enddo
103         write(lunout,*) 'W91 u', MAXVAL(ucov(:,1)), MINVAL(ucov(:,1))
104         ucov(:,1)=ucov(:,1)*cu
105c vitesse meridienne vcov
106         do j=1,jjm
107            sinth  = sin(rlatv(j))
108            costh  = cos(rlatv(j))
109            Ath = -rad*K*R0*sinth*(costh**(R0-1))
110            do i=1,iip1
111               ij=(j-1)*iip1+i
112               lon = rlonv(i)
113               vcov(ij,1) = Ath*sin(R0*lon)
114            enddo
115         enddo
116         write(lunout,*) 'W91 v', MAXVAL(vcov(:,1)), MINVAL(vcov(:,1))
117         vcov(:,1)=vcov(:,1)*cv
118         
119c         ucov=0
120c         vcov=0
121      ELSE
122c test non-tournant, onde se propageant en latitude
123         do j=1,jjp1
124            do i=1,iip1
125               ij=(j-1)*iip1+i
126               ps(ij) = 1e5*(1 + .1*exp(-100*(1+sin(rlatu(j)))**2) )
127            enddo
128         enddo
129         
130c     rho = preff/(cpp*teta)
131         teta = .01*preff/cpp   ! rho = 100 ; phi = ps/rho = 1e3 ; c=30 m/s = 2600 km/j = 23 degres / j
132         ucov=0.
133         vcov=0.
134      END IF     
135     
136      CALL pression ( ip1jmp1, ap, bp, ps, p       )
137      CALL massdair(p,masse)
138
139      END
140c-----------------------------------------------------------------------
Note: See TracBrowser for help on using the repository browser.