psmile_mg_search_3d_real.F90

Go to the documentation of this file.
00001 !-----------------------------------------------------------------------
00002 ! Copyright 2006-2010, NEC Europe Ltd., London, UK.
00003 ! All rights reserved. Use is subject to OASIS4 license terms.
00004 !-----------------------------------------------------------------------
00005 !BOP
00006 !
00007 ! !ROUTINE: PSMILe_mg_search_3d_real
00008 !
00009 ! !INTERFACE:
00010 
00011       subroutine psmile_mg_search_3d_real (comp_info,       &
00012                         found, locations, len, search_data, &
00013                         grid_id, tol, ierror)
00014 !
00015 ! !USES:
00016 !
00017       use PRISM_constants
00018       use psmile_grid, only : common_grid_range
00019       use PSMILe, dummy_interface => PSMILe_mg_search_3d_real
00020 
00021       implicit none
00022 !
00023 ! !INPUT PARAMETERS:
00024 !
00025       Type (Enddef_comp), Intent (In)       :: comp_info
00026 !
00027 !     Info on the component in which the donor cells
00028 !     should be searched.
00029 !
00030       Type (Enddef_search_data)             :: search_data
00031 !
00032 !     Info's on coordinates to be searched
00033 !
00034       Integer, Intent (In)                  :: len (search_data%npart)
00035 !
00036 !     Number of coords to be searched
00037 !
00038       Integer, Intent (In)                  :: grid_id
00039 !
00040 !     Grid id
00041 !
00042       Real, Intent (In)                     :: tol
00043 !
00044 !     Absolute tolerance for search of "identical" points
00045 !     TOL >= 0.0
00046 !
00047 ! !INPUT/OUTPUT PARAMETERS:
00048 !
00049       Type (integer_vector), Intent(InOut)  :: found (search_data%npart)
00050 
00051 !     Finest level number on which a grid cell was found for point I.
00052 !     Level number = -(nlev+1): Never found (input value)
00053 
00054       Type (integer_vector), Intent(InOut)  :: locations (search_data%npart)
00055 !
00056 !     Indices of the grid cell in which the point was found.
00057 !     Assumed input value locations(:)%vector (:, len) = 0
00058 !
00059 ! !OUTPUT PARAMETERS:
00060 !
00061       Integer, Intent (Out)                 :: ierror
00062 
00063 !     Returns the error code of PSMILe_mg_search_3d_real;
00064 !             ierror = 0 : No error
00065 !             ierror > 0 : Severe error
00066 !
00067 ! !LOCAL VARIABLES
00068 !
00069       Type (Corner_Block), Pointer :: corner_pointer
00070       Integer                      :: i, ipart
00071 !
00072 !     ... for locations searched
00073 !
00074       Integer                      :: ibeg, iend
00075       Integer                      :: control (2, ndim_3d, search_data%npart)
00076 !
00077       Type (real_vector)           :: array (ndim_3d, search_data%npart)
00078       Integer                      :: shape (2, ndim_3d, search_data%npart)
00079       Integer                      :: j, k, len1, len2
00080 !
00081 !     ... for levels
00082 !
00083       Integer                      :: lev, levc, nlev
00084       Integer                      :: ijkinc (ndim_3d), ijkcoa (ndim_3d)
00085       Type(Grid), Pointer          :: grid_info
00086 !
00087 !     ... for error parameters
00088 !
00089       Integer, parameter           :: nerrp = 2
00090       Integer                      :: ierrp (nerrp)
00091 !
00092 ! !DESCRIPTION:
00093 !
00094 ! Subroutine "PSMILe_mg_search_3d_real" searches the donor
00095 ! for the subgrid sent by the sending process.
00096 !
00097 ! Subroutine "PSMILe_Search_donor_irreg2_real" is designed for grids of
00098 ! type "PRISM_Irrlonlatvrt". We assume that this routine is only called
00099 ! for methods that are used for coupling.
00100 !
00101 !
00102 ! !REVISION HISTORY:
00103 !
00104 !   Date      Programmer   Description
00105 ! ----------  ----------   -----------
00106 ! 03.07.21    H. Ritzdorf  created
00107 ! 20.04.2011  M. Hanke     moved parts of this routine into
00108 !                          psmile_search_donor_irreg3_dble.F90
00109 !
00110 !EOP
00111 !----------------------------------------------------------------------
00112 !
00113 ! $Id: psmile_mg_search_3d_real.F90 3175 2011-05-04 10:48:06Z hanke $
00114 ! $Author: hanke $
00115 !
00116    Character(len=len_cvs_string), save :: mycvs = 
00117        '$Id: psmile_mg_search_3d_real.F90 3175 2011-05-04 10:48:06Z hanke $'
00118 !
00119 !----------------------------------------------------------------------
00120 !
00121 !  Initialization
00122 !
00123 #ifdef VERBOSE
00124       print 9990, trim(ch_id), comp_info%comp_id
00125 
00126       call psmile_flushstd
00127 #endif /* VERBOSE */
00128 
00129       ierror = 0
00130       grid_info => Grids(grid_id)
00131 !
00132       corner_pointer => Grids(grid_id)%corner_pointer
00133 !
00134 #ifdef PRISM_ASSERTION
00135       if ( grid_info%nlev <= 0 ) then
00136          call psmile_assert (__FILE__, __LINE__, "nlev <= 0")
00137       endif
00138 #endif
00139 !
00140 !-----------------------------------------------------------------------
00141 !     Generate 3d array if necessary
00142 !-----------------------------------------------------------------------
00143 !
00144       if (search_data%grid_type == PRISM_Reglonlatvrt .or. &
00145           search_data%grid_type == PRISM_Irrlonlat_regvrt) then
00146 !        ... Dimensions are only 1-dimensional or 2-dimensional
00147 !        ... Generate 3d arrays
00148 !
00149          shape = search_data%range(:, :, 1:search_data%npart)
00150 !
00151          do ipart = 1, search_data%npart
00152             Allocate (array(1,ipart)%vector(len(ipart)), STAT = ierror)
00153             if ( ierror > 0 ) then
00154                ierrp (1) = ierror
00155                ierrp (2) = len (ipart)
00156                ierror = PRISM_Error_Alloc
00157                call psmile_error ( ierror, 'array (1, ipart)', &
00158                                    ierrp, 2, __FILE__, __LINE__ )
00159                return
00160             endif
00161 !
00162             Allocate (array(2,ipart)%vector(len(ipart)), STAT = ierror)
00163             if ( ierror > 0 ) then
00164                ierrp (1) = ierror
00165             ierrp (2) = len (ipart)
00166                ierror = PRISM_Error_Alloc
00167             call psmile_error ( ierror, 'array2', &
00168                                    ierrp, 2, __FILE__, __LINE__ )
00169                return
00170             endif
00171 !
00172             Allocate (array(3,ipart)%vector(len(ipart)), STAT = ierror)
00173             if ( ierror > 0 ) then
00174                ierrp (1) = ierror
00175             ierrp (2) = len (ipart)
00176                ierror = PRISM_Error_Alloc
00177             call psmile_error ( ierror, 'array3', &
00178                                    ierrp, 2, __FILE__, __LINE__ )
00179                return
00180             endif
00181 !
00182 !===> ... Generate 3-dimensional temporary arrays
00183 !
00184             len1 =  search_data%range(2,1, ipart)-search_data%range(1,1, ipart) + 1
00185             len2 = (search_data%range(2,2, ipart)-search_data%range(1,2, ipart) + 1) * len1
00186 
00187             do k = 1, search_data%range(2,3,ipart)-search_data%range(1,3,ipart) + 1
00188                array(3,ipart)%vector((k-1)*len2+1:k*len2) = &
00189                   search_data%search_real(3, ipart)%vector(k)
00190             end do
00191 
00192             if (search_data%grid_type == PRISM_Reglonlatvrt) then
00193 #ifdef PRISM_ASSERTION
00194                if (search_data%dim_size(1, ipart) /= len1) then
00195                   call psmile_assert (__FILE__, __LINE__, &
00196                                       "dims(1, ipart) != len1")
00197                endif
00198 
00199                if (search_data%dim_size(1, ipart)*search_data%dim_size(2, ipart) /= len2) then
00200                   call psmile_assert (__FILE__, __LINE__, &
00201                                       "dims(1, ipart)*dims(2, ipart) != len2")
00202                endif
00203 #endif
00204 !
00205                do k = 1, search_data%range(2,3, ipart)-search_data%range(1,3, ipart) + 1
00206                   do j = 1, search_data%range(2,2, ipart)-search_data%range(1,2, ipart) + 1
00207                      array(1,ipart)%vector ((k-1)*len2+(j-1)*len1+1:(k-1)*len2+j*len1) = &
00208                         search_data%search_dble(1,ipart)%vector(1:len1)
00209                      array(2,ipart)%vector ((k-1)*len2+(j-1)*len1+1:(k-1)*len2+j*len1) = &
00210                         search_data%search_dble(2,ipart)%vector(j)
00211                   end do
00212                end do
00213 
00214             else if (search_data%grid_type == PRISM_Irrlonlat_regvrt) then
00215 
00216 #ifdef PRISM_ASSERTION
00217                if (search_data%dim_size(1, ipart) /= len2) then
00218                   print *, 'dims(1, ipart)', search_data%dim_size(1, ipart), len2
00219                   call psmile_assert (__FILE__, __LINE__, &
00220                                       "dims(1, ipart) != len2")
00221                endif
00222 
00223                if (search_data%dim_size(2, ipart) /= len2) then
00224                   print *, 'dims(2, ipart)', search_data%dim_size(2, ipart), len2
00225                   call psmile_assert (__FILE__, __LINE__, &
00226                                       "dims(2, ipart) != len2")
00227                endif
00228 #endif
00229 
00230                do k = 1, search_data%range(2,3, ipart)-search_data%range(1,3, ipart) + 1
00231                   do j = 1, len2
00232                      array(1,ipart)%vector ((k-1)*len2+1:k*len2) = &
00233                         search_data%search_real(1, ipart)%vector(1:len2)
00234                      array(2,ipart)%vector ((k-1)*len2+1:k*len2) = &
00235                         search_data%search_real(2, ipart)%vector(1:len2)
00236                   end do
00237                end do
00238             endif
00239 !
00240          end do ! ipart
00241 !
00242       else
00243 !
00244          do ipart = 1, search_data%npart
00245             array(1,ipart)%vector => search_data%search_real(1,ipart)%vector
00246             array(2,ipart)%vector => search_data%search_real(2,ipart)%vector
00247             array(3,ipart)%vector => search_data%search_real(3,ipart)%vector
00248          end do ! ipart
00249 !
00250          shape = search_data%shape(1:2, 1:ndim_3d, 1:search_data%npart)
00251       endif
00252 !
00253 !-----------------------------------------------------------------------
00254 !     Coarsest level nlev
00255 !-----------------------------------------------------------------------
00256 !
00257       control (:, :, 1:search_data%npart) = search_data%range (:, :, 1:search_data%npart)
00258 !
00259       nlev = grid_info%nlev
00260 !
00261       do ipart = 1, search_data%npart
00262          lev = nlev
00263 !
00264          ibeg = 1
00265          iend = len (ipart)
00266 
00267 #ifdef PRISM_ASSERTION
00268          if (grid_info%mg_infos(lev)%levdim(1) /= 0 .or. &
00269              grid_info%mg_infos(lev)%levdim(2) /= 0 .or. &
00270              grid_info%mg_infos(lev)%levdim(3) /= 0) then
00271 
00272             call psmile_assert (__FILE__, __LINE__, &
00273                                 "coarsest level dim != 0")
00274          endif
00275 #endif
00276 !
00277          call psmile_mg_coarse_3d_real (lev, &
00278                      grid_info%mg_infos(lev)%real_arrays%chmin,    &
00279                      grid_info%mg_infos(lev)%real_arrays%chmax,    &
00280                      found(ipart)%vector, locations(ipart)%vector, &
00281                      array(1,ipart)%vector, array(2,ipart)%vector, &
00282                      array(3,ipart)%vector,  &
00283                      ibeg, iend)
00284 !
00285          if (ibeg > iend) then
00286 !           Make control range empty
00287             control (1, 1, ipart) = control (2, 1, ipart)
00288             control (2, 1, ipart) = control (2, 1, ipart) - 1
00289             cycle
00290          endif
00291 !
00292 !===> ... Find range
00293 !  den range besser einschraenken ! 3d Range !
00294 !
00295 !        call psmile_find_range (found(ipart)%vector, ibeg, iend,
00296 !                     control,
00297 !                     iret,  ierror)
00298 !        if (ierror > 0) return
00299 !
00300 !        if (iret == 0) go to 2000
00301 
00302 !
00303 !===> Multiple grids
00304 !
00305 !        ijkinc (1) = max (1, (iend-ibeg)/2)
00306 !        ijkinc (2) = max (1, (jend-jbeg)/2)
00307 !        ijkinc (3) = max (1, (kend-kbeg)/2)
00308 !
00309          ijkinc (:) = 1
00310 !
00311          do lev = nlev-1, 1, -1
00312 !
00313             levc = lev + 1
00314 !
00315             ijkcoa (:) = 2
00316 !
00317             do i = 1, ndim_3d
00318                if (grid_info%mg_infos(lev)%levdim(i) == &
00319                    grid_info%mg_infos(levc)%levdim(i)) then
00320                   ijkcoa (i) = 1
00321                endif
00322             enddo
00323 !
00324 ! chmin, chmax, midp als 4d array alloziieren !
00325 !
00326             call psmile_mg_next_level_3d_real ( grid_id, lev, nlev  ,  &
00327                 grid_info%mg_infos(lev)%real_arrays%chmin(1)%vector,   &
00328                 grid_info%mg_infos(lev)%real_arrays%chmin(2)%vector,   &
00329                 grid_info%mg_infos(lev)%real_arrays%chmin(3)%vector,   &
00330                 grid_info%mg_infos(lev)%real_arrays%chmax(1)%vector,   &
00331                 grid_info%mg_infos(lev)%real_arrays%chmax(2)%vector,   &
00332                 grid_info%mg_infos(lev)%real_arrays%chmax(3)%vector,   &
00333                 grid_info%mg_infos(lev)%real_arrays%midp(1)%vector,    &
00334                 grid_info%mg_infos(lev)%real_arrays%midp(2)%vector,    &
00335                 grid_info%mg_infos(lev)%real_arrays%midp(3)%vector,    &
00336                 grid_info%mg_infos(lev)%levdim,                        &
00337                 found(ipart)%vector, locations(ipart)%vector,          &
00338                 search_data%range(:,:, ipart),                         &
00339                 array(1,ipart)%vector, array(2,ipart)%vector,          &
00340                 array(3,ipart)%vector, shape (:, :, ipart),            &
00341                 control (:, :, ipart), ijkinc, ijkcoa, ierror)
00342 
00343             if (ierror > 0) return
00344 !
00345 !           ijkinc (:) = max (1, ijkinc(:)/ijkcoa(:))
00346          enddo ! lev
00347 !
00348 #ifdef SEARCH_ON_FINAL_GRID
00349 !
00350 !   Get final locations controlling the real cells
00351 !   Glaube nicht, dass das sich lohnt. Zu teuer. Lasse es trotzdem mal drin
00352 !
00353          if (corner_pointer%nbr_corners == 8) then
00354             call psmile_mg_final_level_3d_real ( comp_id, nlev,            &
00355                       found(ipart)%vector,  locations(ipart)%vector,       &
00356                       search_data%range (:, :, ipart),                     &
00357                       array(1,ipart)%vector, array(2,ipart)%vector,        &
00358                       array(3,ipart)%vector, shape (:, :, ipart), control, &
00359                       grid_id,                                             &
00360                       corner_pointer%corners_real(1)%vector,               &
00361                       corner_pointer%corners_real(2)%vector,               &
00362                       corner_pointer%corners_real(3)%vector,               &
00363                       corner_pointer%corner_shape,                         &
00364                       corner_pointer%nbr_corners,                          &
00365                       grid_info%ijk0, tol, ierror)
00366             if (ierror > 0) return
00367          endif
00368 #else
00369 !
00370 !   Compute locations
00371 !
00372 !cdir vector
00373          do i = 1, len (ipart)
00374             locations(ipart)%vector((i-1)*ndim_3d+1) = &
00375             locations(ipart)%vector((i-1)*ndim_3d+1) + grid_info%ijk0 (1)
00376             locations(ipart)%vector((i-1)*ndim_3d+2) = &
00377             locations(ipart)%vector((i-1)*ndim_3d+2) + grid_info%ijk0 (2)
00378             locations(ipart)%vector((i-1)*ndim_3d+3) = &
00379             locations(ipart)%vector((i-1)*ndim_3d+3) + grid_info%ijk0 (3)
00380          end do
00381 !
00382 #endif
00383       end do ! ipart
00384 
00385 !
00386 !-----------------------------------------------------------------------
00387 !     All done; Free memory allocated
00388 !-----------------------------------------------------------------------
00389 !
00390       if (search_data%grid_type == PRISM_Reglonlatvrt .or. &
00391           search_data%grid_type == PRISM_Irrlonlat_regvrt) then
00392 
00393          do ipart = 1, search_data%npart
00394             Deallocate (array(3,ipart)%vector)
00395             Deallocate (array(2,ipart)%vector)
00396             Deallocate (array(1,ipart)%vector)
00397          end do ! ipart
00398       endif
00399 !
00400 #ifdef VERBOSE
00401       print 9980, trim(ch_id), grid_info%comp_id, ierror
00402 
00403       call psmile_flushstd
00404 #endif /* VERBOSE */
00405 !
00406 !  Formats:
00407 !
00408 9990 format (1x, a, ': psmile_mg_search_3d_real: comp_id =', i3)
00409 9980 format (1x, a, ': psmile_mg_search_3d_real: comp_id =', i3, &
00410                     '; eof ierror =', i4)
00411 
00412       end subroutine PSMILe_mg_search_3d_real