[222] | 1 | ! |
---|
| 2 | ! $Header: /home/cvsroot/LMDZ4/libf/phylmd/ocean_slab_mod.F90,v 1.3 2008-02-04 16:24:28 fairhead Exp $ |
---|
| 3 | ! |
---|
| 4 | MODULE surf_heat_transp_mod |
---|
| 5 | |
---|
| 6 | |
---|
| 7 | CONTAINS |
---|
| 8 | |
---|
| 9 | SUBROUTINE divgrad_phy(nlevs,temp,delta) |
---|
| 10 | |
---|
| 11 | |
---|
| 12 | |
---|
| 13 | IMPLICIT NONE |
---|
| 14 | |
---|
| 15 | #include "dimensions.h" |
---|
| 16 | #include "dimphys.h" |
---|
| 17 | #include "paramet.h" |
---|
| 18 | #include "comcstfi.h" |
---|
| 19 | #include "comgeom.h" |
---|
| 20 | #include "comhdiff.h" |
---|
| 21 | |
---|
| 22 | INTEGER nlevs, ll |
---|
| 23 | REAL temp(ngridmx,nlevs) |
---|
| 24 | REAL delta_2d(ip1jmp1,nlevs) |
---|
| 25 | REAL delta(ngridmx,nlevs) |
---|
| 26 | REAL ghx(ip1jmp1,nlevs), ghy(ip1jm,nlevs) |
---|
| 27 | |
---|
| 28 | |
---|
| 29 | CALL gr_fi_dyn(nlevs,ngridmx,iip1,jjp1,temp,delta_2d) |
---|
| 30 | CALL grad(nlevs,delta_2d,ghx,ghy) |
---|
| 31 | DO ll=1,nlevs |
---|
| 32 | ghx(:,ll)=ghx(:,ll)*zmasqu |
---|
| 33 | ! pas de diffusion zonale |
---|
| 34 | ! ghx(:,ll)=0. |
---|
| 35 | ghy(:,ll)=ghy(:,ll)*zmasqv |
---|
| 36 | END DO |
---|
| 37 | CALL diverg(nlevs,ghx,ghy,delta_2d) |
---|
| 38 | CALL gr_dyn_fi(nlevs,iip1,jjp1,ngridmx,delta_2d,delta) |
---|
| 39 | |
---|
| 40 | |
---|
| 41 | END SUBROUTINE divgrad_phy |
---|
| 42 | |
---|
| 43 | |
---|
| 44 | |
---|
| 45 | SUBROUTINE init_masquv(zmasq) |
---|
| 46 | |
---|
| 47 | IMPLICIT NONE |
---|
| 48 | |
---|
| 49 | #include "dimensions.h" |
---|
| 50 | #include "dimphys.h" |
---|
| 51 | #include "paramet.h" |
---|
| 52 | #include "comcstfi.h" |
---|
| 53 | #include "comgeom.h" |
---|
| 54 | #include "comhdiff.h" |
---|
| 55 | |
---|
| 56 | |
---|
| 57 | REAL zmasq(ngridmx) |
---|
| 58 | REAL zmasq_2d(ip1jmp1) |
---|
| 59 | REAL ff(ip1jm) |
---|
| 60 | REAL eps |
---|
| 61 | INTEGER i |
---|
| 62 | |
---|
| 63 | |
---|
| 64 | ! Masques u,v |
---|
| 65 | zmasqu=1. |
---|
| 66 | zmasqv=1. |
---|
| 67 | |
---|
| 68 | CALL gr_fi_dyn(1,ngridmx,iip1,jjp1,zmasq,zmasq_2d) |
---|
| 69 | |
---|
| 70 | DO i=1,ip1jmp1-1 |
---|
| 71 | IF (zmasq_2d(i).GT.1e-5 .OR. zmasq_2d(i+1).GT.1e-5) THEN |
---|
| 72 | zmasqu(i)=0. |
---|
| 73 | ENDIF |
---|
| 74 | END DO |
---|
| 75 | DO i=iip1,ip1jmp1,iip1 |
---|
| 76 | zmasqu(i)=zmasqu(i-iim) |
---|
| 77 | END DO |
---|
| 78 | DO i=1,ip1jm |
---|
| 79 | IF (zmasq_2d(i).GT.1e-5 .OR. zmasq_2d(i+iip1).GT.1e-5) THEN |
---|
| 80 | zmasqv(i)=0. |
---|
| 81 | END IF |
---|
| 82 | END DO |
---|
| 83 | |
---|
| 84 | |
---|
| 85 | ! Coriolis (pour Ekman transp.) |
---|
| 86 | eps=1e-5 |
---|
| 87 | ! CALL getin('slab_eps',eps) |
---|
| 88 | ! print *,'epsilon=',eps |
---|
| 89 | ff=fext*unsairez |
---|
| 90 | DO i=1,ip1jm |
---|
| 91 | unsev(i)=eps/(ff(i)*ff(i)+eps**2) |
---|
| 92 | unsfv(i)=ff(i)/(ff(i)*ff(i)+eps**2) |
---|
| 93 | ENDDO |
---|
| 94 | CALL gr_v_scal(1,unsfv,unsfu) |
---|
| 95 | CALL gr_v_scal(1,unsev,unseu) |
---|
| 96 | ! Alpha variable? |
---|
| 97 | ! alpha_var=.FALSE. |
---|
| 98 | ! CALL getin('slab_alphav',alpha_var) |
---|
| 99 | |
---|
| 100 | |
---|
| 101 | |
---|
| 102 | END SUBROUTINE init_masquv |
---|
| 103 | |
---|
| 104 | |
---|
| 105 | |
---|
| 106 | SUBROUTINE slab_ekman2(tx_phy,ty_phy,ts_phy,dt_phy) |
---|
| 107 | |
---|
| 108 | use slab_ice_h |
---|
| 109 | |
---|
| 110 | IMPLICIT NONE |
---|
| 111 | |
---|
| 112 | #include "dimensions.h" |
---|
| 113 | #include "dimphys.h" |
---|
| 114 | #include "paramet.h" |
---|
| 115 | #include "comcstfi.h" |
---|
| 116 | #include "callkeys.h" |
---|
| 117 | #include "comgeom.h" |
---|
| 118 | #include "comhdiff.h" |
---|
| 119 | |
---|
| 120 | INTEGER ij |
---|
| 121 | REAL txv(ip1jm),fluxm(ip1jm),tyv(ip1jm) |
---|
| 122 | REAL fluxtm(ip1jm,noceanmx),fluxtz(ip1jmp1,noceanmx) |
---|
| 123 | REAL tyu(ip1jmp1),txu(ip1jmp1),fluxz(ip1jmp1),fluxv(ip1jmp1) |
---|
| 124 | REAL dt(ip1jmp1,noceanmx),ts(ip1jmp1,noceanmx) |
---|
| 125 | REAL tx_phy(ngridmx),ty_phy(ngridmx) |
---|
| 126 | REAL dt_phy(ngridmx,noceanmx),ts_phy(ngridmx,noceanmx) |
---|
| 127 | |
---|
| 128 | |
---|
| 129 | |
---|
| 130 | |
---|
| 131 | ! Passage taux,y sur grilles 2d |
---|
| 132 | CALL gr_fi_dyn(1,ngridmx,iip1,jjp1,tx_phy,txu) |
---|
| 133 | CALL gr_fi_dyn(1,ngridmx,iip1,jjp1,ty_phy,tyu) |
---|
| 134 | ! Passage grille u,v |
---|
| 135 | ! Multiplication par f ou eps. |
---|
| 136 | CALL gr_v_scal(1,txu,txv) |
---|
| 137 | CALL gr_v_scal(1,tyu,tyv) |
---|
| 138 | fluxm=tyv*unsev-txv*unsfv |
---|
| 139 | ! fluxm=-txv*unsfv |
---|
| 140 | CALL gr_u_scal(1,txu,txu) |
---|
| 141 | CALL gr_u_scal(1,tyu,tyu) |
---|
| 142 | fluxz=tyu*unsfu+txu*unseu |
---|
| 143 | ! fluxz=tyu*unsfu |
---|
| 144 | |
---|
| 145 | ! Calcul de T, Tdeep |
---|
| 146 | CALL gr_fi_dyn(2,ngridmx,iip1,jjp1,ts_phy,ts) |
---|
| 147 | |
---|
| 148 | ! Flux de masse |
---|
| 149 | fluxm=fluxm*cv*cuvsurcv*zmasqv |
---|
| 150 | fluxz=fluxz*cu*cvusurcu*zmasqu |
---|
| 151 | ! Flux de masse vertical |
---|
| 152 | DO ij=iip2,ip1jm |
---|
| 153 | fluxv(ij)=fluxz(ij)-fluxz(ij-1)-fluxm(ij)+fluxm(ij-iip1) |
---|
| 154 | ENDDO |
---|
| 155 | DO ij=iip1,ip1jmi1,iip1 |
---|
| 156 | fluxv(ij+1)=fluxv(ij+iip1) |
---|
| 157 | END DO |
---|
| 158 | ! Poles |
---|
| 159 | fluxv(1)=-SUM(fluxm(1:iim)) |
---|
| 160 | fluxv(ip1jmp1)=SUM(fluxm(ip1jm-iim:ip1jm-1)) |
---|
| 161 | fluxv=fluxv |
---|
| 162 | |
---|
| 163 | !Calcul flux de chaleur méridiens |
---|
| 164 | DO ij=1,ip1jm |
---|
| 165 | fluxtm(ij,1)=fluxm(ij)*(ts(ij+iip1,1)+ts(ij,1))/2. |
---|
| 166 | fluxtm(ij,2)=-fluxm(ij)*(ts(ij+iip1,2)+ts(ij,2))/2. |
---|
| 167 | END DO |
---|
| 168 | |
---|
| 169 | !Calcul flux chaleur zonaux |
---|
| 170 | DO ij=iip2,ip1jm |
---|
| 171 | IF (fluxz(ij).GE.0.) THEN |
---|
| 172 | fluxtz(ij,1)=fluxz(ij)*ts(ij,1) |
---|
| 173 | fluxtz(ij,2)=-fluxz(ij)*ts(ij+1,2) |
---|
| 174 | ELSE |
---|
| 175 | fluxtz(ij,1)=fluxz(ij)*ts(ij+1,1) |
---|
| 176 | fluxtz(ij,2)=-fluxz(ij)*ts(ij,2) |
---|
| 177 | ENDIF |
---|
| 178 | END DO |
---|
| 179 | DO ij=iip1*2,ip1jmp1,iip1 |
---|
| 180 | fluxtz(ij,:)=fluxtz(ij-iim,:) |
---|
| 181 | END DO |
---|
| 182 | |
---|
| 183 | ! Calcul de dT |
---|
| 184 | DO ij=iip2,ip1jm |
---|
| 185 | dt(ij,:)=fluxtz(ij-1,:)-fluxtz(ij,:) & |
---|
| 186 | +fluxtm(ij,:)-fluxtm(ij-iip1,:) |
---|
| 187 | IF (fluxv(ij).GT.0.) THEN |
---|
| 188 | dt(ij,1)=dt(ij,1)+fluxv(ij)*ts(ij,2) |
---|
| 189 | dt(ij,2)=dt(ij,2)-fluxv(ij)*ts(ij,2) |
---|
| 190 | ELSE |
---|
| 191 | dt(ij,1)=dt(ij,1)+fluxv(ij)*ts(ij,1) |
---|
| 192 | dt(ij,2)=dt(ij,2)-fluxv(ij)*ts(ij,1) |
---|
| 193 | ENDIF |
---|
| 194 | dt(ij,:)=dt(ij,:)/aire(ij) |
---|
| 195 | END DO |
---|
| 196 | DO ij=iip1,ip1jmi1,iip1 |
---|
| 197 | dt(ij+1,:)=dt(ij+iip1,:) |
---|
| 198 | END DO |
---|
| 199 | ! Pôles |
---|
| 200 | dt(1,:)=SUM(fluxtm(1:iim,:),dim=1) |
---|
| 201 | IF (fluxv(1).GT.0.) THEN |
---|
| 202 | dt(1,1)=dt(1,1)+fluxv(1)*ts(1,2) |
---|
| 203 | dt(1,2)=dt(1,2)-fluxv(1)*ts(1,2) |
---|
| 204 | ELSE |
---|
| 205 | dt(1,1)=dt(1,1)+fluxv(1)*ts(1,1) |
---|
| 206 | dt(1,2)=dt(1,2)-fluxv(1)*ts(1,1) |
---|
| 207 | ENDIF |
---|
| 208 | dt(1,:)=dt(1,:)/apoln |
---|
| 209 | dt(ip1jmp1,:)=-SUM(fluxtm(ip1jm-iim:ip1jm-1,:),dim=1) |
---|
| 210 | IF (fluxv(ip1jmp1).GT.0.) THEN |
---|
| 211 | dt(ip1jmp1,1)=dt(ip1jmp1,1)+fluxv(ip1jmp1)*ts(ip1jmp1,2) |
---|
| 212 | dt(ip1jmp1,2)=dt(ip1jmp1,2)-fluxv(ip1jmp1)*ts(ip1jmp1,2) |
---|
| 213 | ELSE |
---|
| 214 | dt(ip1jmp1,1)=dt(ip1jmp1,1)+fluxv(ip1jmp1)*ts(ip1jmp1,1) |
---|
| 215 | dt(ip1jmp1,2)=dt(ip1jmp1,2)-fluxv(ip1jmp1)*ts(ip1jmp1,1) |
---|
| 216 | ENDIF |
---|
| 217 | dt(ip1jmp1,:)=dt(ip1jmp1,:)/apols |
---|
| 218 | dt(2:iip1,1)=dt(1,1) |
---|
| 219 | dt(2:iip1,2)=dt(1,2) |
---|
| 220 | dt(ip1jm+1:ip1jmp1,1)=dt(ip1jmp1,1) |
---|
| 221 | dt(ip1jm+1:ip1jmp1,2)=dt(ip1jmp1,2) |
---|
| 222 | |
---|
| 223 | ! Retour grille physique |
---|
| 224 | CALL gr_dyn_fi(2,iip1,jjp1,ngridmx,dt,dt_phy) |
---|
| 225 | |
---|
| 226 | |
---|
| 227 | END SUBROUTINE slab_ekman2 |
---|
| 228 | |
---|
| 229 | |
---|
| 230 | END MODULE surf_heat_transp_mod |
---|
| 231 | |
---|
| 232 | |
---|
| 233 | |
---|