psmile_def_domains.F90

Go to the documentation of this file.
00001 !-----------------------------------------------------------------------
00002 ! Copyright 2007-2010, SGI Germany, Munich, Germany.
00003 ! All rights reserved. Use is subject to OASIS4 license terms.
00004 !-----------------------------------------------------------------------
00005 !BOP
00006 !
00007 ! !ROUTINE: PSMILE_Def_Domains
00008 !
00009 ! !INTERFACE:
00010 !
00011       subroutine psmile_def_domains ( ierror )
00012 !
00013 ! !USES:
00014 !
00015       use PSMILe, dummy_interface => PSMILE_Def_Domains
00016       implicit none
00017 
00018       include 'prism.inc'
00019 !
00020 ! !OUTPUT PARAMETERS:
00021 !
00022       integer, Intent (Out)               :: ierror
00023 
00024 !
00025 ! !LOCAL VARIABLES:
00026 !
00027       integer                 :: ii,kk,jj,ll,iiloc,il_i
00028       integer                 :: il_taskid
00029       integer                             :: ierrp(2)
00030       integer                             :: iinbr
00031       integer                             :: max_nbr
00032       integer                 :: iicomp_id
00033       integer                 :: iicomp_comm
00034       integer                 :: iicomp_size
00035       integer                             :: iirecvcount
00036       integer                             :: iisendcount
00037       integer                 :: iigrid_id,iigrid_type
00038       integer                             :: len,iindim
00039       integer,pointer                     :: iioffsets(:,:,:)
00040       integer,pointer                     :: iiextblks(:,:,:)
00041       integer,pointer                     :: iioffsetx(:)
00042       integer,pointer                     :: iioffsety(:)
00043       integer,pointer                     :: iishapes(:,:,:)
00044       integer,pointer                     :: iitmp_shape(:)
00045       integer,pointer                     :: iiextents(:,:)
00046       integer,pointer                     :: nnbr(:)
00047       integer                             :: iiglob_extents(3)
00048       integer                             :: iilayout(2)
00049       integer                             :: minijk(3)
00050       logical                             :: ll_overlap
00051       logical                             :: is_regular
00052       logical,pointer                     :: ll_mask(:,:)
00053       integer                 :: noc,nvertices(3),nof(3),lenz
00054 !
00055 !
00056 ! !DESCRIPTION:
00057 !
00058 !  Driver for mpp_def_domains. Loops of fields which are associated
00059 !  with file I/O and computes the global extent of a grid. 
00060 !  The extends of the task local partitions and their offsets in global
00061 !  index space are fed into mpp_def_domains. A check is performed on the
00062 !  grid's valid shape if indices are passed as taks local indices.
00063 !  If the test is true the user has to provide additional  informations
00064 !  regarding the partitioning via prism_def_partition.
00065 !
00066 ! !REVISION HISTORY:
00067 !
00068 !   Date      Programmer   Description
00069 ! ----------  ----------   -----------
00070 ! 9.12.03    R. Vogelsang created
00071 ! 23.3.04    R. Vogelsang Added the code to handle  multiple I/O tasks
00072 ! 02.10.05   R. Vogelsang Added extensions for writing of vertical bounds 
00073 !             and handling of Gaussian reduced grids.
00074 ! 28.02.06   S. Valcke    Copyright string modified
00075 ! 1.03.06    R. Vogelsang Final synchronization with local repository
00076 !
00077 !EOP
00078 ! ---------------------------------------------------------------------
00079 
00080 Character(len=len_cvs_string), save :: mycvs = 
00081        '$Id: psmile_def_domains.F90 2687 2010-10-28 15:15:52Z coquart $'
00082 
00083       ierror = 0
00084 #ifdef __PSMILE_WITH_IO
00085 #ifdef VERBOSE
00086       print*,trim(ch_id),' : psmile_def_domains: start'
00087       call psmile_flushstd
00088 #endif
00089 
00090       do ii=1,Number_of_fields_allocated
00091 !RV   If true we have a file subject to I/O
00092         if(associated(Fields(ii)%io_chan_infos)) then
00093 !
00094 !       Start loop over I/O channels, means transients subject to file I/O!
00095 !
00096         do il_taskid=1,size(Fields(ii)%io_task_lookup)
00097 
00098           il_i=Fields(ii)%io_task_lookup(il_taskid)
00099 
00100           if (il_i.gt.0) then
00101 !
00102 !           Let the I/O infos point ot the actual I/O channel infos
00103 !
00104             Fields(ii)%io_infos =>  Fields(ii)%io_chan_infos(il_i)
00105 
00106           else
00107             CYCLE
00108           endif
00109           if(il_i .eq. 1 ) then
00110 
00111 !
00112           iicomp_id=Fields(ii)%comp_id
00113 !rv          iicomp_comm=Comps(iicomp_id)%comm
00114           iicomp_comm=Comps(iicomp_id)%act_comm
00115 !rv          iicomp_size=Comps(iicomp_id)%size
00116           call MPI_Comm_size(iicomp_comm,iicomp_size,ierror)
00117 
00118           if(.not.associated(IO_Comps_infos(iicomp_id)%pelist)) then
00119              ierror=PRISM_Error_IO_Domain
00120              ierrp(1)=iicomp_id
00121              call psmile_error ( ierror &
00122                               ,'pelist for component does not exist!' &
00123                               ,ierrp, 1, __FILE__,__LINE__)
00124           endif
00125 #ifdef VERBOSE
00126           print*,trim(ch_id),  &
00127           ' : psmile_def_domains: The pelist of component',iicomp_id,' is ', &
00128           IO_Comps_infos(iicomp_id)%pelist
00129           call psmile_flushstd
00130 #endif
00131           call mpp_set_current_pelist(IO_Comps_infos(iicomp_id)%pelist)
00132 !
00133           iigrid_id=Methods(Fields(ii)%method_id)%grid_id
00134 !RV       For bounds 
00135           noc=4
00136           if(associated(Grids(iigrid_id)%corner_pointer)) &
00137              noc=Grids(iigrid_id)%corner_pointer%nbr_corners
00138 #ifdef VERBOSE
00139           print*,trim(ch_id), & 
00140           ' : psmile_def_domains: Number of corners: ',noc
00141           call psmile_flushstd
00142 #endif
00143 
00144 !RV       For bounds 
00145           iigrid_type=Grids(iigrid_id)%grid_type
00146           if ( iigrid_type == PRISM_unstructlonlatvrt) then
00147              len = 1
00148           else if ( iigrid_type == PRISM_unstructlonlat_sigmavrt .or. &
00149                     iigrid_type == PRISM_unstructlonlat_regvrt   ) then
00150              len = 2
00151           else
00152              len = ndim_3d
00153           endif
00154 !RV
00155 !RV       For vertical bounds
00156 !RV
00157         iindim=Grids(iigrid_id)%n_dim
00158         lenz=Grids(iigrid_id)%grid_shape(2,iindim) - &
00159              Grids(iigrid_id)%grid_shape(1,iindim) + 1
00160 
00161         is_regular=.true.
00162         select case ( Grids(iigrid_id)%grid_type )
00163 
00164         case ( PRISM_Irrlonlat_sigmavrt)
00165 
00166            nof(1:2) = noc/2
00167            nof(ndim_3d) = 2
00168 
00169            do jj=1,3
00170              nvertices(jj)=noc
00171            enddo
00172            nvertices(3)=nof(3)
00173 
00174         case ( PRISM_Irrlonlat_regvrt )
00175 
00176            nof(1:2)     = noc/2
00177            nof(ndim_3d) = 2
00178 
00179            do jj=1,3
00180              nvertices(jj)=noc
00181            enddo
00182            nvertices(3)=nof(3)
00183 
00184         case ( PRISM_Reglonlat_sigmavrt )
00185 
00186            nof(1:ndim_3d) = 2
00187            do jj=1,3
00188              nvertices(jj)=noc
00189            enddo
00190            nvertices(3)=nof(3)
00191 
00192         case ( PRISM_Reglonlatvrt )
00193 
00194            nof(1:ndim_3d) = 2
00195 
00196            do jj=1,3
00197              nvertices(jj)=noc
00198            enddo
00199            nvertices(3)=nof(3)
00200 
00201         case ( PRISM_Irrlonlatvrt )
00202 
00203            nof(1:ndim_3d) = noc
00204 
00205            do jj=1,3
00206              nvertices(jj)=noc
00207            enddo
00208         case ( PRISM_Unstructlonlatvrt )
00209 
00210            nof(1:ndim_3d) = noc
00211 
00212            do jj=1,3
00213              nvertices(jj)=noc
00214            enddo
00215            is_regular=.false.
00216 
00217         case ( PRISM_Unstructlonlat_regvrt )
00218 
00219            nof(1:2) = noc/2
00220            nof(ndim_3d) = 2
00221 
00222            do jj=1,3
00223              nvertices(jj)=noc
00224            enddo
00225 
00226            nvertices(3)=nof(3)
00227            is_regular=.false.
00228 
00229            ierrp (1) = Grids(iigrid_id)%grid_type
00230 
00231            ierror = PRISM_Error_Internal
00232 
00233            call psmile_error(ierror,'grid generation type not supported yet', &
00234                              ierrp, 1, __FILE__, __LINE__ )
00235            return
00236 
00237         case ( PRISM_Unstructlonlat_sigmavrt )
00238 
00239            nof(1:2) = noc/2
00240            nof(ndim_3d) = 2
00241 
00242            do jj=1,3
00243              nvertices(jj)=noc
00244            enddo
00245            nvertices(3)=nof(3)
00246            is_regular=.false.
00247 
00248            ierrp (1) = Grids(iigrid_id)%grid_type
00249 
00250            ierror = PRISM_Error_Internal
00251 
00252            call psmile_error(ierror,'grid generation type not supported yet', &
00253                              ierrp, 1, __FILE__, __LINE__ )
00254            return
00255 
00256         case ( PRISM_Gaussreduced_regvrt )
00257   
00258 !rv to be discussed with Rene
00259            nof(1:2) = noc/2
00260            nof(ndim_3d) = 2
00261 
00262            do jj=1,3
00263              nvertices(jj)=noc
00264            enddo
00265   
00266            nvertices(3)=nof(3)
00267            is_regular=.false.
00268 
00269   
00270         case ( PRISM_Gaussreduced_sigmavrt )
00271   
00272            nof(1:2) = noc/2
00273            nof(ndim_3d) = 2
00274 
00275            do jj=1,3
00276              nvertices(jj)=noc
00277            enddo
00278            nvertices(3)=nof(3)
00279            is_regular=.false.
00280   
00281         case DEFAULT
00282 
00283            ierrp (1) = Grids(iigrid_id)%grid_type
00284 
00285            ierror = PRISM_Error_Internal
00286 
00287            call psmile_error(ierror, 'unknown grid generation type', &
00288                                ierrp, 1, __FILE__, __LINE__ )
00289            return
00290 
00291         end select
00292 
00293 
00294 
00295 !RV          if (.not.associated(Grids( iigrid_id)%partition)) then
00296 !RV             ierror=PRISM_Error_IO_Domain
00297 !RV             ierrp(1)=iigrid_id
00298 !RV             ierrp(2)=iigrid_type
00299 !RV             call PSMILe_Error ( ierror &
00300 !RV                           ,'For I/O prism_def_partition has to be called!' &
00301 !RV                              ,ierrp, 2, __FILE__,__LINE__)
00302 !RV          endif
00303 
00304 
00305           allocate(Fields(ii)%io_infos%p_mpp_io,stat=ierror)
00306           if (ierror > 0) then
00307              ierrp (1) = ierror
00308              ierrp (2) = 1
00309              ierror = PRISM_Error_Alloc
00310              call psmile_error ( ierror, 'p_mpp_io', &
00311                              ierrp, 2, __FILE__, __LINE__ )
00312              return
00313           endif
00314 !
00315 !         Nullify of substructures and pointers
00316 !
00317           Nullify(Fields(ii)%io_infos%p_mpp_io%domain)
00318           Nullify(Fields(ii)%io_infos%p_mpp_io%xdom)
00319           Nullify(Fields(ii)%io_infos%p_mpp_io%ydom)
00320           Nullify(Fields(ii)%io_infos%p_mpp_io%zdom)
00321 !          
00322           Nullify(Fields(ii)%io_infos%p_mpp_io%ig_vars)
00323           Nullify(Fields(ii)%io_infos%p_mpp_io%field)
00324           Nullify(Fields(ii)%io_infos%p_mpp_io%latij)
00325           Nullify(Fields(ii)%io_infos%p_mpp_io%lonij)
00326           Nullify(Fields(ii)%io_infos%p_mpp_io%vertij)
00327           Nullify(Fields(ii)%io_infos%p_mpp_io%maskij)
00328           Nullify(Fields(ii)%io_infos%p_mpp_io%crnij)
00329           Nullify(Fields(ii)%io_infos%p_mpp_io%subaij)
00330 !
00331           Nullify(Fields(ii)%io_infos%p_mpp_io%x_axis)
00332           Nullify(Fields(ii)%io_infos%p_mpp_io%y_axis)
00333           Nullify(Fields(ii)%io_infos%p_mpp_io%z_axis)
00334           Nullify(Fields(ii)%io_infos%p_mpp_io%t_axis)
00335           Nullify(Fields(ii)%io_infos%p_mpp_io%b_axis)
00336           Nullify(Fields(ii)%io_infos%p_mpp_io%blk_axis)
00337           Nullify(Fields(ii)%io_infos%p_mpp_io%anaxis)
00338           Nullify(Fields(ii)%io_infos%p_mpp_io%crnaxis)
00339 !
00340           Nullify(Fields(ii)%io_infos%p_mpp_io%atts)
00341 !
00342           allocate(Fields(ii)%io_infos%p_mpp_io%domain(4),stat=ierror)
00343           if (ierror > 0) then
00344              ierrp (1) = ierror
00345              ierrp (2) = 1
00346              ierror = PRISM_Error_Alloc
00347              call psmile_error ( ierror, 'domain', &
00348                              ierrp, 2, __FILE__, __LINE__ )
00349              return
00350           endif
00351           allocate(Fields(ii)%io_infos%p_mpp_io%xdom(4),stat=ierror)
00352           if (ierror > 0) then
00353              ierrp (1) = ierror
00354              ierrp (2) = 1
00355              ierror = PRISM_Error_Alloc
00356              call psmile_error ( ierror, 'xdom', &
00357                              ierrp, 2, __FILE__, __LINE__ )
00358              return
00359           endif
00360           allocate(Fields(ii)%io_infos%p_mpp_io%ydom(4),stat=ierror)
00361           if (ierror > 0) then
00362              ierrp (1) = ierror
00363              ierrp (2) = 1
00364              ierror = PRISM_Error_Alloc
00365              call psmile_error ( ierror, 'ydom', &
00366                              ierrp, 2, __FILE__, __LINE__ )
00367              return
00368           endif
00369           allocate(Fields(ii)%io_infos%p_mpp_io%zdom(4),stat=ierror)
00370           if (ierror > 0) then
00371              ierrp (1) = ierror
00372              ierrp (2) = 1
00373              ierror = PRISM_Error_Alloc
00374              call psmile_error ( ierror, 'zdom', &
00375                              ierrp, 2, __FILE__, __LINE__ )
00376              return
00377           endif
00378 
00379           allocate(iishapes(2,1:Grids(iigrid_id)%n_dim,1:iicomp_size) &
00380                            ,stat=ierror)
00381 
00382           if (ierror > 0) then
00383              ierrp (1) = ierror
00384              ierrp (2) = 1
00385              ierror = PRISM_Error_Alloc
00386              call psmile_error ( ierror, 'iishapes', &
00387                              ierrp, 2, __FILE__, __LINE__ )
00388              return
00389           endif
00390 
00391           allocate(iiextents(1:Grids(iigrid_id)%n_dim,1:iicomp_size) &
00392                   ,stat=ierror)
00393 
00394           if (ierror > 0) then
00395              ierrp (1) = ierror
00396              ierrp (2) = 1
00397              ierror = PRISM_Error_Alloc
00398              call psmile_error ( ierror, 'iextents', &
00399                              ierrp, 2, __FILE__, __LINE__ )
00400              return
00401           endif
00402 
00403 !         iinbr should be always 1. Anyway,find iinbr for compatability
00404 !         reasons.
00405           iinbr=1
00406           max_nbr=iinbr
00407           allocate(nnbr(1:iicomp_size),stat=ierror)
00408 
00409           if (ierror > 0) then
00410              ierrp (1) = ierror
00411              ierrp (2) = 1
00412              ierror = PRISM_Error_Alloc
00413              call psmile_error ( ierror, 'nnbr', &
00414                              ierrp, 2, __FILE__, __LINE__ )
00415              return
00416           endif
00417 
00418           nnbr(1:iicomp_size)=iinbr
00419 
00420           if(associated(Grids( iigrid_id)%partition)) then
00421             iinbr=size(Grids( iigrid_id)%partition,DIM=1)
00422             call MPI_Allgather(iinbr &
00423                             ,1,MPI_INTEGER &
00424                             ,nnbr,1,MPI_INTEGER &
00425                             ,iicomp_comm, ierror)
00426 
00427             do il_i=1,iicomp_size
00428               max_nbr=max(nnbr(il_i),max_nbr)
00429             enddo
00430           endif
00431 #ifdef VERBOSE
00432       print*,trim(ch_id),' : psmile_def_domains: Max. no. of blocks ',max_nbr
00433       call psmile_flushstd
00434 #endif
00435         
00436 
00437           allocate(iioffsets(max_nbr,1:len,1:iicomp_size),stat=ierror)
00438 
00439           if (ierror > 0) then
00440              ierrp (1) = ierror
00441              ierrp (2) = 1
00442              ierror = PRISM_Error_Alloc
00443              call psmile_error ( ierror, 'iioffsets', &
00444                              ierrp, 2, __FILE__, __LINE__ )
00445              return
00446           endif
00447 
00448           allocate(iiextblks(max_nbr,1:len,1:iicomp_size),stat=ierror)
00449           if (ierror > 0) then
00450              ierrp (1) = ierror
00451              ierrp (2) = 1
00452              ierror = PRISM_Error_Alloc
00453              call psmile_error ( ierror, 'iiextblks', &
00454                              ierrp, 2, __FILE__, __LINE__ )
00455              return
00456           endif
00457 
00458           allocate(iioffsetx(1:iicomp_size),stat=ierror)
00459 
00460           if (ierror > 0) then
00461              ierrp (1) = ierror
00462              ierrp (2) = 1
00463              ierror = PRISM_Error_Alloc
00464              call psmile_error ( ierror, 'iioffsetx', &
00465                              ierrp, 2, __FILE__, __LINE__ )
00466              return
00467           endif
00468 
00469           allocate(iioffsety(1:iicomp_size),stat=ierror)
00470 
00471           if (ierror > 0) then
00472              ierrp (1) = ierror
00473              ierrp (2) = 1
00474              ierror = PRISM_Error_Alloc
00475              call psmile_error ( ierror, 'iioffsety', &
00476                              ierrp, 2, __FILE__, __LINE__ )
00477              return
00478           endif
00479 
00480           allocate(ll_mask(iicomp_size,iicomp_size),stat=ierror)
00481 
00482           if (ierror > 0) then
00483              ierrp (1) = ierror
00484              ierrp (2) = 1
00485              ierror = PRISM_Error_Alloc
00486              call psmile_error ( ierror, 'll_mask', &
00487                              ierrp, 2, __FILE__, __LINE__ )
00488              return
00489           endif
00490           
00491 !RV   Collect shapes and global offsets.
00492          
00493           iindim=Grids(iigrid_id)%n_dim
00494           iirecvcount=iindim*2
00495           iisendcount=iindim*2
00496           allocate(iitmp_shape(iirecvcount),stat=ierror)
00497           if (ierror > 0) then
00498              ierrp (1) = ierror
00499              ierrp (2) = 1
00500              ierror = PRISM_Error_Alloc
00501              call psmile_error ( ierror, 'iitmp_shape', &
00502                              ierrp, 2, __FILE__, __LINE__ )
00503              return
00504           endif
00505 
00506           iiloc=0
00507           do kk=1,iindim
00508             do jj=1,2
00509               iiloc=iiloc+1
00510               iitmp_shape(iiloc)=Grids(iigrid_id)%grid_shape(jj,kk)  
00511             enddo
00512           enddo
00513 
00514           call MPI_Allgather(iitmp_shape &
00515                             ,iisendcount,MPI_INTEGER &
00516                             ,iishapes(1,1,1),iirecvcount,MPI_INTEGER &
00517                             ,iicomp_comm, ierror)
00518 
00519           deallocate(iitmp_shape)
00520 
00521 !
00522 !         The user has called prism_def_partition.
00523 !         Offsets declared thre are going to be used. Collect them!
00524 !
00525           if(associated(Grids(iigrid_id)%partition)) then
00526 #ifdef VERBOSE
00527       print*,trim(ch_id),' : psmile_def_domains: Collecting block offsets'
00528       call psmile_flushstd
00529 #endif
00530 
00531           iirecvcount=max_nbr*len
00532           iisendcount=max_nbr*len
00533 
00534           allocate(iitmp_shape(iirecvcount),stat=ierror)
00535           if (ierror > 0) then
00536              ierrp (1) = ierror
00537              ierrp (2) = 1
00538              ierror = PRISM_Error_Alloc
00539              call psmile_error ( ierror, 'iitmp_shape', &
00540                              ierrp, 2, __FILE__, __LINE__ )
00541              return
00542           endif
00543 
00544           iiloc=0
00545           iitmp_shape(:)=0
00546           do kk=1,len
00547             do jj=1,iinbr
00548               iiloc=iiloc+1
00549               iitmp_shape(iiloc)=Grids(iigrid_id)%partition(jj,kk)
00550             enddo
00551           enddo
00552          
00553           call MPI_Allgather(iitmp_shape &
00554                             ,iisendcount,MPI_INTEGER &
00555                             ,iioffsets(1,1,1),iirecvcount,MPI_INTEGER &
00556                             ,iicomp_comm, ierror)
00557           if(iinbr.gt.1) then
00558           iiloc=0
00559           iitmp_shape(:)=0
00560           do kk=1,len
00561             do jj=1,iinbr
00562               iiloc=iiloc+1
00563               iitmp_shape(iiloc)=Grids(iigrid_id)%extent(jj,kk)
00564             enddo
00565           enddo
00566           call MPI_Allgather(iitmp_shape &
00567                             ,iisendcount,MPI_INTEGER &
00568                             ,iiextblks(1,1,1),iirecvcount,MPI_INTEGER &
00569                             ,iicomp_comm, ierror)
00570           endif
00571 
00572           deallocate(iitmp_shape)
00573 !
00574           else
00575 !
00576 !         The user has  N O T  called prism_def_partition.
00577 !         I assume that the lower indices of the shapes are going to give 
00578 !         me the proper  offsets .
00579 !
00580 #ifdef VERBOSE
00581       print*,trim(ch_id),' : psmile_def_domains: Computing block offsets'
00582       call psmile_flushstd
00583 #endif
00584 !         Perform a check if the blocks overlap.  
00585 !         If that is the case the user has declared potentially
00586 !         the shapes in a local index space. He must call prism_def_partition.
00587 !
00588           do kk=1,iicomp_size-1
00589             do jj=kk+1,iicomp_size
00590               ll_overlap=.true.
00591               do ll=1,len
00592               if(.not.(((iishapes(1,ll,jj).ge.iishapes(1,ll,kk)) .and. &
00593                  (iishapes(1,ll,jj).le.iishapes(2,ll,kk))) .or.  &
00594                  ((iishapes(2,ll,jj).ge.iishapes(1,ll,kk)) .and. &
00595                  (iishapes(2,ll,jj).le.iishapes(2,ll,kk))) .and. &
00596                   ll_overlap)) then
00597                   ll_overlap=.false.
00598               endif
00599               enddo
00600               if(ll_overlap) then
00601                 ierror=PRISM_Error_IO_Domain
00602                 ierrp(1)=iigrid_id
00603                 ierrp(2)=iigrid_type
00604                 call psmile_error ( ierror &
00605                            ,'For I/O prism_def_partition has to be called!' &
00606                            ,ierrp, 2, __FILE__,__LINE__)
00607 
00608               endif
00609             enddo
00610           enddo
00611 !
00612 !         Compute offsets
00613 !
00614           minijk(1:3)=2**30
00615           do kk=1,len
00616             do jj=1,iicomp_size
00617               minijk(kk)=min(minijk(kk),iishapes(1,kk,jj))
00618             enddo
00619           enddo
00620 !
00621           do kk=1,len
00622             do jj=1,iicomp_size
00623               iioffsets(1,kk,jj)=iishapes(1,kk,jj)-minijk(kk)
00624             enddo
00625           enddo
00626           endif
00627 !
00628 !RV       Compute extents of all blocks
00629 !
00630           if(iinbr.eq.1) then
00631           do jj=1,iindim
00632             do kk=1,iicomp_size
00633               iiextents(jj,kk)=iishapes(2,jj,kk)-iishapes(1,jj,kk)+1
00634             enddo
00635 #ifdef VERBOSE
00636       print*,trim(ch_id) &
00637      ,' : psmile_def_domains: comp id,varid,grid id' &
00638      ,',direction : valid extents(1:<size of component>)' &
00639      ,iicomp_id,',',ii,',',iigrid_id,',',jj,':',iiextents(jj,1:iicomp_size)
00640       call psmile_flushstd
00641 #endif 
00642           enddo
00643           else
00644           if( iigrid_type == PRISM_Gaussreduced_regvrt) then
00645             do kk=1,iicomp_size
00646               iiextents(1,kk)=sum(iiextblks(1:nnbr(kk),1,kk))
00647               iiextents(2,kk)=1
00648               iiextents(3,kk)=iiextblks(1,3,kk)
00649             enddo
00650 #ifdef VERBOSE
00651             do jj=1,iindim
00652                print*,trim(ch_id) &
00653                     ,' : psmile_def_domains: comp id,varid,grid id' &
00654                     ,',direction : valid extents(1:<size of component>)' &
00655                     ,iicomp_id,',',nnbr,',',iigrid_id,',',jj,':',iiextents(jj,1:iicomp_size)
00656                !      do kk = 1, iicomp_size
00657                !      print*,'The blocks are :',iiextblks(1:nnbr(kk),jj,kk)
00658                !      enddo
00659                call psmile_flushstd
00660             enddo
00661 #endif 
00662           else
00663             do jj=1,iindim
00664               do kk=1,iicomp_size
00665                 iiextents(jj,kk)=sum(iiextblks(1:nnbr(kk),jj,kk))
00666               enddo
00667 #ifdef VERBOSE
00668               print*,trim(ch_id) &
00669                    ,' : psmile_def_domains: comp id,varid,grid id' &
00670                    ,',direction : valid extents(1:<size of component>)' &
00671                    ,iicomp_id,',',nnbr,',',iigrid_id,',',jj,':',iiextents(jj,1:iicomp_size)
00672                !      do kk = 1, iicomp_size
00673                !      print*,'The blocks are :',iiextblks(1:nnbr(kk),jj,kk)
00674                !      enddo
00675               call psmile_flushstd
00676 #endif 
00677             enddo
00678           endif
00679 !RV
00680 !RV       I am overwriting the offsets provided by prism_def_partition!
00681 !RV       Handling of different block offsets per block not yet possible.
00682 !RV
00683           do jj=1,iindim
00684             iioffsets(1,jj,1)=0
00685             do kk=2,iicomp_size
00686               iioffsets(1,jj,kk)=iioffsets(1,jj,kk-1) + iiextents(jj,kk-1)
00687             enddo
00688           enddo
00689           endif
00690 
00691           if(len .eq. 1) then
00692 !RV       Safety dance: If len=1 grid is irregular in all directions.
00693 !RV                     All extents except the first one are 1
00694               iiextents(2,1:iicomp_size)=1
00695               iiextents(3,1:iicomp_size)=1
00696               iioffsets(1,2,1:iicomp_size)=0
00697               iioffsets(1,3,1:iicomp_size)=0
00698           else if (len.eq.2) then
00699 !RV       Safety dance: If len=2 grid is irregular in the horizontal directions.
00700 !RV                     Extents for 2nd direction is one.
00701               iiextents(2,1:iicomp_size)=1
00702               iioffsets(1,2,1:iicomp_size)=0
00703           endif
00704 
00705 !RV       Compute the global extent of the  grid
00706           iiglob_extents=1
00707           if( iigrid_type == PRISM_Gaussreduced_regvrt) then
00708             do kk=1,iicomp_size
00709                iiglob_extents(1)=max(iiglob_extents(1) &
00710                                  ,iioffsets(1,1,kk)+iiextents(1,kk))
00711                iiglob_extents(2)=1
00712                iiglob_extents(3)=max(iiglob_extents(3) &
00713                                  ,iioffsets(1,3,kk)+iiextents(3,kk))
00714             enddo
00715 #ifdef VERBOSE
00716       print*,trim(ch_id) &
00717      ,' : psmile_def_domains: no. of dim,glob. extent of the grid(no. of dim)' &
00718      ,len,',',iiglob_extents(1:3)
00719       call psmile_flushstd
00720 #endif 
00721           else
00722           do jj=1,len
00723              do kk=1,iicomp_size
00724                iiglob_extents(jj)=max(iiglob_extents(jj) &
00725                                  ,iioffsets(1,jj,kk)+iiextents(jj,kk))
00726              enddo
00727           enddo
00728 #ifdef VERBOSE
00729       print*,trim(ch_id) &
00730      ,' : psmile_def_domains: no. of dim,glob. extent of the grid(no. of dim)' &
00731      ,len,',',iiglob_extents(1:len)
00732       call psmile_flushstd
00733 #endif 
00734           endif
00735 
00736           if(len.le.2) then
00737             iioffsetx(1:iicomp_size)=iioffsets(1,1,1:iicomp_size)
00738 !            iioffsety(1:iicomp_size)=(/(jj,jj=0,iicomp_size-1)/)
00739              iioffsety(1:iicomp_size)=0
00740 !            iiglob_extents(2)=iicomp_size
00741           else
00742             iioffsetx(1:iicomp_size)=iioffsets(1,1,1:iicomp_size)
00743             iioffsety(1:iicomp_size)=iioffsets(1,2,1:iicomp_size)
00744           endif
00745 #ifdef VERBOSE
00746       print*,trim(ch_id) &
00747      ,' : psmile_def_domains: offsets in x',iioffsetx(1:iicomp_size)
00748       print*,trim(ch_id) &
00749      ,' : psmile_def_domains: offsets in y',iioffsety(1:iicomp_size)
00750       call psmile_flushstd
00751 #endif 
00752 
00753 !RV       Computation of the layout of the processor grid is hopeless
00754 !RV       for a potentially block structured grid. I am using the block
00755 !RV       approach right from the beginning.
00756 !RV       WARNING: This scheme is not implmented in the MPP_IO package
00757 !RV       distributed with MOM4.
00758           ll_mask=.false.
00759           if(is_regular) then
00760             do jj=1,iicomp_size
00761               ll_mask(jj,jj)=.true.
00762             enddo
00763           else
00764             do jj=1,iicomp_size
00765               ll_mask(jj,1)=.true.
00766             enddo
00767           endif
00768           iilayout(1)=iicomp_size
00769           iilayout(2)=iicomp_size
00770           call mpp_define_domains((/1,iiglob_extents(1),1,iiglob_extents(2)/) &
00771                                  ,iilayout &
00772                                  ,Fields(ii)%io_infos%p_mpp_io%domain(1) &
00773                                  ,xextent=iiextents(1,1:iicomp_size) &
00774                                  ,yextent=iiextents(2,1:iicomp_size) &
00775                                  ,pelist=IO_Comps_infos(iicomp_id)%pelist &
00776                                  ,maskmap=ll_mask &
00777                                  ,offsetx=iioffsetx &
00778                                  ,offsety=iioffsety)
00779 !
00780 !RV       For bounds
00781 !         Lon and Lat
00782           call mpp_define_domains( &
00783                                  (/1,noc*iiglob_extents(1),1,iiglob_extents(2)/) &
00784                                  ,iilayout &
00785                                  ,Fields(ii)%io_infos%p_mpp_io%domain(2) &
00786                              ,xextent=(/(noc*iiextents(1,jj),jj=1,iicomp_size)/) &
00787                              ,yextent=iiextents(2,1:iicomp_size) &
00788                              ,pelist=IO_Comps_infos(iicomp_id)%pelist &
00789                              ,maskmap=ll_mask &
00790                              ,offsetx=(/(noc*iioffsetx(jj),jj=1,iicomp_size)/) &
00791                              ,offsety=iioffsety)
00792 !
00793 !RV       Vertical
00794 !
00795 !          call mpp_define_domains( &
00796 !                          (/1,noc*lenz*iiglob_extents(1),1,iiglob_extents(2)/) &
00797 !                          ,iilayout &
00798 !                          ,Fields(ii)%io_infos%p_mpp_io%domain(3) &
00799 !                          ,xextent=(/(nof(iindim)*lenz*iiextents(1,jj),jj=1,iicomp_size)/) &
00800 !                          ,yextent=iiextents(2,1:iicomp_size) &
00801 !                          ,pelist=IO_Comps_infos(iicomp_id)%pelist &
00802 !                          ,maskmap=ll_mask &
00803 !                        ,offsetx=(/(noc*lenz*iioffsetx(jj),jj=1,iicomp_size)/) &
00804 !                        ,offsety=iioffsety)
00805           call mpp_define_domains( &
00806                           (/1,nof(iindim)*lenz*iiglob_extents(1),1,iiglob_extents(2)/) &
00807                           ,iilayout &
00808                           ,Fields(ii)%io_infos%p_mpp_io%domain(3) &
00809                           ,xextent=(/(nof(iindim)*lenz*iiextents(1,jj),jj=1,iicomp_size)/) &
00810                           ,yextent=iiextents(2,1:iicomp_size) &
00811                           ,pelist=IO_Comps_infos(iicomp_id)%pelist &
00812                           ,maskmap=ll_mask &
00813                         ,offsetx=(/(nof(iindim)*lenz*iioffsetx(jj),jj=1,iicomp_size)/) &
00814                         ,offsety=iioffsety)
00815 !
00816 !RV       For bounds
00817 !
00818           call mpp_get_domain_components(Fields(ii)%io_infos%p_mpp_io%domain(1)&
00819                                         ,Fields(ii)%io_infos%p_mpp_io%xdom(1) &
00820                                         ,Fields(ii)%io_infos%p_mpp_io%ydom(1))
00821 !         Cleanup
00822           deallocate(iishapes,stat=ierror)
00823           if (ierror > 0) then
00824              ierrp (1) = ierror
00825              ierrp (2) = 1
00826              ierror = PRISM_Error_Alloc
00827              call psmile_error ( ierror, 'iishapes', &
00828                              ierrp, 2, __FILE__, __LINE__ )
00829              return
00830           endif
00831 
00832           deallocate(iiextents,stat=ierror)
00833           if (ierror > 0) then
00834              ierrp (1) = ierror
00835              ierrp (2) = 1
00836              ierror = PRISM_Error_Alloc
00837              call psmile_error ( ierror, 'iextents', &
00838                              ierrp, 2, __FILE__, __LINE__ )
00839              return
00840           endif
00841 
00842           deallocate(iioffsets,stat=ierror)
00843           if (ierror > 0) then
00844              ierrp (1) = ierror
00845              ierrp (2) = 1
00846              ierror = PRISM_Error_Alloc
00847              call psmile_error ( ierror, 'iioffsets', &
00848                              ierrp, 2, __FILE__, __LINE__ )
00849              return
00850           endif
00851 
00852           deallocate(iioffsetx,stat=ierror)
00853           if (ierror > 0) then
00854              ierrp (1) = ierror
00855              ierrp (2) = 1
00856              ierror = PRISM_Error_Alloc
00857              call psmile_error ( ierror, 'iioffsetx', &
00858                              ierrp, 2, __FILE__, __LINE__ )
00859              return
00860           endif
00861 
00862           deallocate(iioffsety,stat=ierror)
00863           if (ierror > 0) then
00864              ierrp (1) = ierror
00865              ierrp (2) = 1
00866              ierror = PRISM_Error_Alloc
00867              call psmile_error ( ierror, 'iioffsety', &
00868                              ierrp, 2, __FILE__, __LINE__ )
00869              return
00870           endif
00871 
00872           deallocate(ll_mask,stat=ierror)
00873           if (ierror > 0) then
00874              ierrp (1) = ierror
00875              ierrp (2) = 1
00876              ierror = PRISM_Error_Alloc
00877              call psmile_error ( ierror, 'll_mask', &
00878                              ierrp, 2, __FILE__, __LINE__ )
00879              return
00880           endif
00881 
00882           else  ! For higher IO channels use the domain infos of the first.
00883 
00884            Fields(ii)%io_chan_infos(il_i)%p_mpp_io  => &
00885                     Fields(ii)%io_chan_infos(1)%p_mpp_io
00886           endif
00887 
00888         enddo   ! End of loop over tasks
00889         endif
00890 
00891       enddo
00892 #ifdef VERBOSE
00893       print*,trim(ch_id),' : psmile_def_domains: end'
00894       call psmile_flushstd
00895 #endif
00896 
00897 #endif
00898 
00899       end subroutine PSMILe_Def_Domains