00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018 module ice_restoring
00019
00020
00021
00022 use ice_kinds_mod
00023 use ice_blocks
00024 use ice_calendar, only: dt
00025 use ice_domain
00026 use ice_domain_size
00027 use ice_communicate, only: my_task, master_task
00028 use ice_constants
00029 use ice_exit
00030 use ice_fileunits
00031 use ice_forcing, only: trestore, trest
00032 use ice_state
00033 use ice_timers
00034
00035
00036
00037 implicit none
00038 save
00039
00040 logical (kind=log_kind) ::
00041 restore_ice
00042
00043
00044
00045
00046
00047 real (kind=dbl_kind), dimension (:,:,:,:), allocatable ::
00048 aicen_rest ,
00049 vicen_rest ,
00050 vsnon_rest ,
00051 eicen_rest ,
00052 esnon_rest
00053
00054 real (kind=dbl_kind), dimension (:,:,:,:,:), allocatable ::
00055 trcrn_rest
00056
00057
00058
00059 contains
00060
00061
00062
00063
00064
00065
00066 subroutine ice_HaloRestore_init
00067
00068
00069
00070
00071
00072
00073
00074
00075
00076
00077 allocate (aicen_rest(nx_block,ny_block,ncat,max_blocks), &
00078 vicen_rest(nx_block,ny_block,ncat,max_blocks), &
00079 vsnon_rest(nx_block,ny_block,ncat,max_blocks), &
00080 eicen_rest(nx_block,ny_block,ntilyr,max_blocks), &
00081 esnon_rest(nx_block,ny_block,ntslyr,max_blocks), &
00082 trcrn_rest(nx_block,ny_block,ntrcr,ncat,max_blocks))
00083
00084
00085
00086
00087 call ice_timer_start(timer_bound)
00088 call bound_state (aicen, trcrn, &
00089 vicen, vsnon, &
00090 eicen, esnon)
00091 call ice_timer_stop(timer_bound)
00092
00093 aicen_rest(:,:,:,:) = aicen(:,:,:,:)
00094 vicen_rest(:,:,:,:) = vicen(:,:,:,:)
00095 vsnon_rest(:,:,:,:) = vsnon(:,:,:,:)
00096 eicen_rest(:,:,:,:) = eicen(:,:,:,:)
00097 esnon_rest(:,:,:,:) = esnon(:,:,:,:)
00098 trcrn_rest(:,:,:,:,:) = trcrn(:,:,1:ntrcr,:,:)
00099
00100 if (my_task == master_task) &
00101 write (nu_diag,*) 'ice restoring timescale = ',trestore,' days'
00102
00103 end subroutine ice_HaloRestore_init
00104
00105
00106
00107
00108
00109
00110
00111 subroutine ice_HaloRestore
00112
00113
00114
00115
00116
00117
00118
00119
00120
00121
00122
00123
00124 use ice_distribution
00125
00126
00127
00128
00129
00130
00131
00132
00133
00134 integer (int_kind) ::
00135 i,j,iblk,nt,n,
00136 ilo,ihi,jlo,jhi,
00137 ibc,
00138 npad
00139
00140 type (block) ::
00141 this_block
00142
00143 real (dbl_kind) ::
00144 ctime
00145
00146 call ice_timer_start(timer_bound)
00147
00148
00149
00150
00151
00152
00153
00154
00155 if (trestore == 0) then
00156 trest = dt
00157 else
00158 trest = real(trestore,kind=dbl_kind) * secday
00159 endif
00160 ctime = dt/trest
00161
00162
00163
00164
00165
00166
00167
00168 do iblk = 1, nblocks
00169 this_block = get_block(blocks_ice(iblk),iblk)
00170 ilo = this_block%ilo
00171 ihi = this_block%ihi
00172 jlo = this_block%jlo
00173 jhi = this_block%jhi
00174
00175 if (this_block%iblock == 1) then
00176 if (trim(ew_boundary_type) /= 'cyclic') then
00177 do n = 1, ncat
00178 do j = 1, ny_block
00179 do i = 1, ilo
00180 aicen(i,j,n,iblk) = aicen(i,j,n,iblk) &
00181 + (aicen_rest(i,j,n,iblk)-aicen(i,j,n,iblk))*ctime
00182 vicen(i,j,n,iblk) = vicen(i,j,n,iblk) &
00183 + (vicen_rest(i,j,n,iblk)-vicen(i,j,n,iblk))*ctime
00184 vsnon(i,j,n,iblk) = vsnon(i,j,n,iblk) &
00185 + (vsnon_rest(i,j,n,iblk)-vsnon(i,j,n,iblk))*ctime
00186 do nt = 1, ntrcr
00187 trcrn(i,j,nt,n,iblk) = trcrn(i,j,nt,n,iblk) &
00188 + (trcrn_rest(i,j,nt,n,iblk)-trcrn(i,j,nt,n,iblk))*ctime
00189 enddo
00190 enddo
00191 enddo
00192 enddo
00193 do n = 1, ntilyr
00194 do j = 1, ny_block
00195 do i = 1, ilo
00196 eicen(i,j,n,iblk) = eicen(i,j,n,iblk) &
00197 + (eicen_rest(i,j,n,iblk)-eicen(i,j,n,iblk))*ctime
00198 enddo
00199 enddo
00200 enddo
00201 do n = 1, ntslyr
00202 do j = 1, ny_block
00203 do i = 1, ilo
00204 esnon(i,j,n,iblk) = esnon(i,j,n,iblk) &
00205 + (esnon_rest(i,j,n,iblk)-esnon(i,j,n,iblk))*ctime
00206 enddo
00207 enddo
00208 enddo
00209 endif
00210
00211 elseif (this_block%iblock == nblocks_x) then
00212 if (trim(ew_boundary_type) /= 'cyclic') then
00213
00214 ibc = nx_block + 1
00215 npad = 0
00216 do i = nx_block, 1, - 1
00217 if (this_block%i_glob(i) == 0) then
00218 do j = 1, ny_block
00219 npad = npad + this_block%j_glob(j)
00220 enddo
00221 endif
00222 if (npad == 0) ibc = ibc - 1
00223 enddo
00224
00225 do n = 1, ncat
00226 do j = 1, ny_block
00227 do i = ihi, ibc
00228 aicen(i,j,n,iblk) = aicen(i,j,n,iblk) &
00229 + (aicen_rest(i,j,n,iblk)-aicen(i,j,n,iblk))*ctime
00230 vicen(i,j,n,iblk) = vicen(i,j,n,iblk) &
00231 + (vicen_rest(i,j,n,iblk)-vicen(i,j,n,iblk))*ctime
00232 vsnon(i,j,n,iblk) = vsnon(i,j,n,iblk) &
00233 + (vsnon_rest(i,j,n,iblk)-vsnon(i,j,n,iblk))*ctime
00234 do nt = 1, ntrcr
00235 trcrn(i,j,nt,n,iblk) = trcrn(i,j,nt,n,iblk) &
00236 + (trcrn_rest(i,j,nt,n,iblk)-trcrn(i,j,nt,n,iblk))*ctime
00237 enddo
00238 enddo
00239 enddo
00240 enddo
00241 do n = 1, ntilyr
00242 do j = 1, ny_block
00243 do i = ihi, ibc
00244 eicen(i,j,n,iblk) = eicen(i,j,n,iblk) &
00245 + (eicen_rest(i,j,n,iblk)-eicen(i,j,n,iblk))*ctime
00246 enddo
00247 enddo
00248 enddo
00249 do n = 1, ntslyr
00250 do j = 1, ny_block
00251 do i = ihi, ibc
00252 esnon(i,j,n,iblk) = esnon(i,j,n,iblk) &
00253 + (esnon_rest(i,j,n,iblk)-esnon(i,j,n,iblk))*ctime
00254 enddo
00255 enddo
00256 enddo
00257 endif
00258 endif
00259
00260 if (this_block%jblock == 1) then
00261 if (trim(ns_boundary_type) /= 'cyclic') then
00262 do n = 1, ncat
00263 do j = 1, jlo
00264 do i = 1, nx_block
00265 aicen(i,j,n,iblk) = aicen(i,j,n,iblk) &
00266 + (aicen_rest(i,j,n,iblk)-aicen(i,j,n,iblk))*ctime
00267 vicen(i,j,n,iblk) = vicen(i,j,n,iblk) &
00268 + (vicen_rest(i,j,n,iblk)-vicen(i,j,n,iblk))*ctime
00269 vsnon(i,j,n,iblk) = vsnon(i,j,n,iblk) &
00270 + (vsnon_rest(i,j,n,iblk)-vsnon(i,j,n,iblk))*ctime
00271 do nt = 1, ntrcr
00272 trcrn(i,j,nt,n,iblk) = trcrn(i,j,nt,n,iblk) &
00273 + (trcrn_rest(i,j,nt,n,iblk)-trcrn(i,j,nt,n,iblk))*ctime
00274 enddo
00275 enddo
00276 enddo
00277 enddo
00278 do n = 1, ntilyr
00279 do j = 1, jlo
00280 do i = 1, nx_block
00281 eicen(i,j,n,iblk) = eicen(i,j,n,iblk) &
00282 + (eicen_rest(i,j,n,iblk)-eicen(i,j,n,iblk))*ctime
00283 enddo
00284 enddo
00285 enddo
00286 do n = 1, ntslyr
00287 do j = 1, jlo
00288 do i = 1, nx_block
00289 esnon(i,j,n,iblk) = esnon(i,j,n,iblk) &
00290 + (esnon_rest(i,j,n,iblk)-esnon(i,j,n,iblk))*ctime
00291 enddo
00292 enddo
00293 enddo
00294 endif
00295
00296 elseif (this_block%jblock == nblocks_y) then
00297 if (trim(ns_boundary_type) /= 'cyclic' .and. &
00298 trim(ns_boundary_type) /= 'tripole' .and. &
00299 trim(ns_boundary_type) /= 'tripoleT') then
00300
00301 ibc = ny_block + 1
00302 npad = 0
00303 do j = ny_block, 1, - 1
00304 if (this_block%j_glob(j) == 0) then
00305 do i = 1, nx_block
00306 npad = npad + this_block%i_glob(i)
00307 enddo
00308 endif
00309 if (npad == 0) ibc = ibc - 1
00310 enddo
00311
00312 do n = 1, ncat
00313 do j = jhi, ibc
00314 do i = 1, nx_block
00315 aicen(i,j,n,iblk) = aicen(i,j,n,iblk) &
00316 + (aicen_rest(i,j,n,iblk)-aicen(i,j,n,iblk))*ctime
00317 vicen(i,j,n,iblk) = vicen(i,j,n,iblk) &
00318 + (vicen_rest(i,j,n,iblk)-vicen(i,j,n,iblk))*ctime
00319 vsnon(i,j,n,iblk) = vsnon(i,j,n,iblk) &
00320 + (vsnon_rest(i,j,n,iblk)-vsnon(i,j,n,iblk))*ctime
00321 do nt = 1, ntrcr
00322 trcrn(i,j,nt,n,iblk) = trcrn(i,j,nt,n,iblk) &
00323 + (trcrn_rest(i,j,nt,n,iblk)-trcrn(i,j,nt,n,iblk))*ctime
00324 enddo
00325 enddo
00326 enddo
00327 enddo
00328 do n = 1, ntilyr
00329 do j = jhi, ibc
00330 do i = 1, nx_block
00331 eicen(i,j,n,iblk) = eicen(i,j,n,iblk) &
00332 + (eicen_rest(i,j,n,iblk)-eicen(i,j,n,iblk))*ctime
00333 enddo
00334 enddo
00335 enddo
00336 do n = 1, ntslyr
00337 do j = jhi, ibc
00338 do i = 1, nx_block
00339 esnon(i,j,n,iblk) = esnon(i,j,n,iblk) &
00340 + (esnon_rest(i,j,n,iblk)-esnon(i,j,n,iblk))*ctime
00341 enddo
00342 enddo
00343 enddo
00344 endif
00345 endif
00346
00347 enddo
00348
00349 call ice_timer_stop(timer_bound)
00350
00351 end subroutine ice_HaloRestore
00352
00353
00354
00355 end module ice_restoring
00356
00357