00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026 module ice_ocean
00027
00028
00029
00030 use ice_kinds_mod
00031 use ice_constants
00032
00033
00034
00035 implicit none
00036 save
00037
00038 logical (kind=log_kind) ::
00039 oceanmixed_ice
00040
00041 real (kind=dbl_kind), parameter ::
00042 cprho = cp_ocn*rhow
00043
00044
00045
00046 contains
00047
00048
00049
00050
00051
00052
00053
00054
00055
00056
00057
00058
00059
00060
00061
00062
00063
00064 subroutine ocean_mixed_layer (dt)
00065
00066
00067
00068 use ice_blocks
00069 use ice_domain
00070 use ice_state, only: aice, uvel, vvel
00071 use ice_flux
00072 use ice_grid, only: tmask
00073 use ice_atmo, only: atmo_boundary_layer, atmbndy, atmo_boundary_const
00074
00075
00076
00077 real (kind=dbl_kind), intent(in) ::
00078 dt
00079
00080
00081
00082 real (kind=dbl_kind) ::
00083 TsfK ,
00084 swabs
00085
00086 real (kind=dbl_kind), parameter ::
00087 frzmlt_max = c1000
00088
00089 integer (kind=int_kind) ::
00090 i, j ,
00091 ij ,
00092 iblk ,
00093 ilo,ihi,jlo,jhi
00094
00095 real (kind=dbl_kind), dimension(nx_block,ny_block) ::
00096 delt ,
00097 delq ,
00098 shcoef,
00099 lhcoef
00100
00101 integer (kind=int_kind), save ::
00102 icells
00103
00104 integer (kind=int_kind), dimension(nx_block*ny_block), save ::
00105 indxi, indxj
00106
00107 type (block) ::
00108 this_block
00109
00110 do iblk = 1, nblocks
00111
00112
00113
00114
00115
00116
00117 icells = 0
00118 do j = 1, ny_block
00119 do i = 1, nx_block
00120 if (tmask(i,j,iblk)) then
00121 icells = icells + 1
00122 indxi(icells) = i
00123 indxj(icells) = j
00124 else
00125 sst (i,j,iblk) = c0
00126 frzmlt (i,j,iblk) = c0
00127 flwout_ocn(i,j,iblk) = c0
00128 fsens_ocn (i,j,iblk) = c0
00129 flat_ocn (i,j,iblk) = c0
00130 evap_ocn (i,j,iblk) = c0
00131 endif
00132 enddo
00133 enddo
00134
00135
00136
00137
00138
00139 if (trim(atmbndy) == 'constant') then
00140 call atmo_boundary_const (nx_block, ny_block, &
00141 'ice', icells, &
00142 indxi, indxj, &
00143 uatm (:,:,iblk), &
00144 vatm (:,:,iblk), &
00145 wind (:,:,iblk), &
00146 rhoa (:,:,iblk), &
00147 strairx_ocn(:,:,iblk), &
00148 strairy_ocn(:,:,iblk), &
00149 lhcoef (:,:), &
00150 shcoef (:,:) )
00151 else
00152 call atmo_boundary_layer (nx_block, ny_block, &
00153 'ocn', icells, &
00154 indxi, indxj, &
00155 sst (:,:,iblk), &
00156 potT (:,:,iblk), &
00157 uatm (:,:,iblk), &
00158 vatm (:,:,iblk), &
00159 uvel (:,:,iblk), &
00160 vvel (:,:,iblk), &
00161 wind (:,:,iblk), &
00162 zlvl (:,:,iblk), &
00163 Qa (:,:,iblk), &
00164 rhoa (:,:,iblk), &
00165 strairx_ocn(:,:,iblk), &
00166 strairy_ocn(:,:,iblk), &
00167 Tref_ocn (:,:,iblk), &
00168 Qref_ocn (:,:,iblk), &
00169 delt (:,:), &
00170 delq (:,:), &
00171 lhcoef (:,:), &
00172 shcoef (:,:) )
00173 endif
00174
00175
00176
00177
00178
00179
00180 alvdr_ocn(:,:,iblk) = albocn
00181 alidr_ocn(:,:,iblk) = albocn
00182 alvdf_ocn(:,:,iblk) = albocn
00183 alidf_ocn(:,:,iblk) = albocn
00184
00185
00186
00187
00188
00189
00190
00191
00192 do ij = 1, icells
00193 i = indxi(ij)
00194 j = indxj(ij)
00195
00196
00197 swabs = (c1-alvdr_ocn(i,j,iblk)) * swvdr(i,j,iblk) &
00198 + (c1-alidr_ocn(i,j,iblk)) * swidr(i,j,iblk) &
00199 + (c1-alvdf_ocn(i,j,iblk)) * swvdf(i,j,iblk) &
00200 + (c1-alidf_ocn(i,j,iblk)) * swidf(i,j,iblk)
00201
00202
00203 TsfK = sst(i,j,iblk) + Tffresh
00204
00205
00206 flwout_ocn(i,j,iblk) = -stefan_boltzmann * TsfK**4
00207
00208
00209 fsens_ocn(i,j,iblk) = shcoef(i,j) * delt(i,j)
00210 flat_ocn (i,j,iblk) = lhcoef(i,j) * delq(i,j)
00211 evap_ocn (i,j,iblk) = -flat_ocn(i,j,iblk) / Lvap
00212
00213
00214 sst(i,j,iblk) = sst(i,j,iblk) + dt * ( &
00215 (fsens_ocn(i,j,iblk) + flat_ocn(i,j,iblk) + flwout_ocn(i,j,iblk) &
00216 + flw(i,j,iblk) + swabs) * (c1-aice(i,j,iblk)) &
00217 + fhocn(i,j,iblk) + fswthru(i,j,iblk)) &
00218 / (cprho*hmix(i,j,iblk))
00219
00220
00221 if (sst(i,j,iblk) <= Tf(i,j,iblk) .and. qdp(i,j,iblk) > c0) qdp(i,j,iblk) = c0
00222
00223
00224 sst(i,j,iblk) = sst(i,j,iblk) - qdp(i,j,iblk)*dt/(cprho*hmix(i,j,iblk))
00225
00226
00227 frzmlt(i,j,iblk) = (Tf(i,j,iblk)-sst(i,j,iblk))*cprho*hmix(i,j,iblk)/dt
00228 frzmlt(i,j,iblk) = min(max(frzmlt(i,j,iblk),-frzmlt_max),frzmlt_max)
00229
00230
00231 if (sst(i,j,iblk) <= Tf(i,j,iblk)) sst(i,j,iblk) = Tf(i,j,iblk)
00232
00233 enddo
00234 enddo
00235
00236 end subroutine ocean_mixed_layer
00237
00238
00239
00240 end module ice_ocean
00241
00242