[12] | 1 | MODULE etat0_williamson_mod |
---|
[19] | 2 | USE icosa |
---|
[12] | 3 | PRIVATE |
---|
[204] | 4 | REAL(rstd), PARAMETER :: h0=8.E3 |
---|
| 5 | REAL(rstd), PARAMETER :: R0=4 |
---|
| 6 | REAL(rstd), PARAMETER :: K0=7.848E-6 |
---|
| 7 | REAL(rstd), PARAMETER :: omega0=K0 |
---|
[12] | 8 | |
---|
[204] | 9 | PUBLIC getin_etat0, compute_etat0 |
---|
[12] | 10 | |
---|
[204] | 11 | CONTAINS |
---|
[12] | 12 | |
---|
[204] | 13 | SUBROUTINE getin_etat0 |
---|
| 14 | USE mpipara, ONLY : is_mpi_root |
---|
| 15 | USE disvert_mod, ONLY : caldyn_eta, eta_lag |
---|
[164] | 16 | IF(caldyn_eta /= eta_lag) THEN |
---|
[204] | 17 | IF(is_mpi_root) PRINT *, 'etat0_type=williamson91.6 (Williamson,1991) must be used with caldyn_eta=eta_lag' |
---|
[164] | 18 | STOP |
---|
| 19 | END IF |
---|
| 20 | |
---|
| 21 | IF(llm>1) THEN |
---|
[204] | 22 | IF(is_mpi_root) PRINT *, 'etat0_type=williamson91.6 (Williamson,1991) must be used with llm=1 but llm =',llm |
---|
[164] | 23 | STOP |
---|
| 24 | END IF |
---|
[204] | 25 | END SUBROUTINE getin_etat0 |
---|
[164] | 26 | |
---|
[204] | 27 | SUBROUTINE compute_etat0(ngrid,lon,lat, phis,mass,rhodz,ulon,ulat) |
---|
| 28 | IMPLICIT NONE |
---|
| 29 | INTEGER, INTENT(IN) :: ngrid |
---|
| 30 | REAL(rstd),INTENT(IN) :: lon(ngrid) |
---|
| 31 | REAL(rstd),INTENT(IN) :: lat(ngrid) |
---|
| 32 | REAL(rstd),INTENT(OUT) :: phis(ngrid) |
---|
| 33 | REAL(rstd),INTENT(OUT) :: mass(ngrid) |
---|
| 34 | REAL(rstd),INTENT(OUT) :: rhodz(ngrid) |
---|
| 35 | REAL(rstd),INTENT(OUT) :: ulon(ngrid) |
---|
| 36 | REAL(rstd),INTENT(OUT) :: ulat(ngrid) |
---|
[164] | 37 | |
---|
[204] | 38 | REAL(rstd) :: coslat,sinlat, A,B,C |
---|
| 39 | INTEGER :: ij |
---|
[12] | 40 | |
---|
[204] | 41 | DO ij=1,ngrid |
---|
| 42 | coslat=cos(lat(ij)) |
---|
| 43 | sinlat=sin(lat(ij)) |
---|
| 44 | A = 0.5*omega0*(2*omega+omega0)*coslat**2 & |
---|
| 45 | + 0.25*K0**2*coslat**(2*R0)*((R0+1)*coslat**2+(2*R0**2-R0-2)-2*R0**2/coslat**2) |
---|
| 46 | B = 2*(omega+omega0)*K0/((R0+1)*(R0+2))*coslat**R0*((R0**2+2*R0+2)-(R0+1)**2*coslat**2) |
---|
| 47 | C = 0.25*K0**2*coslat**(2*R0)*((R0+1)*coslat**2-(R0+2)) |
---|
[12] | 48 | |
---|
[204] | 49 | mass(ij) = (g*h0+radius**2*A+radius**2*B*cos(R0*lon(ij))+radius**2*C*cos(2*R0*lon(ij)))/g |
---|
| 50 | ulon(ij) = radius*omega0*coslat+radius*K0*coslat**(R0-1)*(R0*sinlat**2-coslat**2)*cos(R0*lon(ij)) |
---|
| 51 | ulat(ij) = -radius*K0*R0*coslat**(R0-1)*sinlat*sin(R0*lon(ij)) |
---|
[12] | 52 | ENDDO |
---|
| 53 | |
---|
[204] | 54 | phis(:)=0. |
---|
| 55 | rhodz(:)=mass(:) |
---|
[12] | 56 | |
---|
[204] | 57 | END SUBROUTINE compute_etat0 |
---|
[12] | 58 | |
---|
| 59 | END MODULE etat0_williamson_mod |
---|