00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011 subroutine psmile_ccompact_gauss2_dble ( send_info, &
00012 grid_valid_shape, shape, nb_corners, &
00013 array_x, array_y, array_z, &
00014 extra_search, dest_size, nbr_cells_tot, &
00015 source_cell_index, &
00016 neighcells, dest_x, dest_y, dest_z, ierror )
00017
00018
00019
00020 use PRISM_constants
00021
00022 use PSMILe, dummy_interface => PSMILe_ccompact_gauss2_dble
00023
00024 implicit none
00025
00026
00027
00028 Type (Send_information), Intent (InOut) :: send_info
00029
00030
00031
00032 Integer, Intent (In) :: grid_valid_shape (2, ndim_3d)
00033
00034
00035
00036 Integer, Intent (In) :: shape (2, ndim_3d)
00037
00038
00039
00040 Integer, Intent (In) :: nb_corners
00041
00042
00043
00044 Double Precision, Intent (In) :: array_x ( shape(1,1):shape(2,1),
00045 nb_corners )
00046 Double Precision, Intent (In) :: array_y ( shape(1,1):shape(2,1),
00047 nb_corners )
00048 Double Precision, Intent (In) :: array_z ( shape(1,3):shape(2,3) )
00049
00050
00051
00052 Type (Extra_search_info), Intent (InOut) :: extra_search
00053
00054
00055
00056 Integer, Intent (In) :: dest_size
00057
00058
00059
00060 Integer, Intent (In) :: nbr_cells_tot
00061
00062
00063
00064 Integer, Intent (InOut) :: source_cell_index (nbr_cells_tot)
00065
00066
00067
00068 Integer, Intent (InOut) :: neighcells (nbr_cells_tot, nb_corners)
00069
00070
00071
00072
00073
00074 Double Precision, Intent (Out) :: dest_x (2*dest_size)
00075 Double Precision, Intent (Out) :: dest_y (2*dest_size)
00076 Double Precision, Intent (Out) :: dest_z (2*dest_size)
00077
00078
00079
00080 Integer, Intent (Out) :: ierror
00081
00082
00083
00084
00085
00086
00087
00088
00089
00090 Integer, Pointer :: list_entries (:, :)
00091 Integer :: nbr_cells_loc
00092
00093
00094
00095 Integer :: ni, nij, nij_loc
00096 Integer :: i, k, l, n
00097
00098
00099
00100 Integer :: dest_idx_off
00101 Integer :: stride, idx1, idx2
00102
00103
00104
00105
00106
00107
00108
00109
00110
00111
00112
00113
00114
00115
00116
00117
00118
00119
00120
00121
00122
00123
00124
00125
00126
00127
00128 Character(len=len_cvs_string), save :: mycvs =
00129 '$Id: psmile_ccompact_gauss2_dble.F90 2590 2010-09-23 15:56:33Z hanke $'
00130
00131
00132
00133
00134
00135 #ifdef VERBOSE
00136 print 9990, trim(ch_id), send_info%send_entire_valid_shape
00137
00138 call psmile_flushstd
00139 #endif /* VERBOSE */
00140
00141 #ifdef PRISM_ASSERTION
00142 if (dest_size < send_info%n_list-send_info%num2recv) then
00143 print *, trim(ch_id), 'dest_size, n_list, num2recv', &
00144 dest_size, send_info%n_list, send_info%num2recv
00145 call psmile_assert (__FILE__, __LINE__, &
00146 "dest_size is not sufficient")
00147 endif
00148 #endif
00149
00150 ierror = 0
00151 nbr_cells_loc = nbr_cells_tot - send_info%num2recv
00152
00153
00154
00155
00156
00157
00158
00159 if (send_info%send_entire_valid_shape) then
00160
00161 ni = grid_valid_shape (2, 1) - grid_valid_shape (1, 1) + 1
00162 nij = 2*ni
00163
00164 do i = grid_valid_shape (1, 1), grid_valid_shape (2, 1)
00165 dest_x (i) = min(array_x(i,1),array_x(i,2))
00166 dest_y (i) = min(array_y(i,1),array_y(i,2))
00167 end do
00168
00169 do i = grid_valid_shape (1, 1), grid_valid_shape (2, 1)
00170 dest_x (ni+i) = max(array_x(i,1),array_x(i,2))
00171 dest_y (ni+i) = max(array_y(i,1),array_y(i,2))
00172 end do
00173
00174 dest_z (1:nij) = array_z(grid_valid_shape(1,3))
00175
00176
00177
00178
00179 do k = 2, grid_valid_shape (2, 3) - grid_valid_shape (1, 3) + 1
00180 dest_x((k-1)*nij+1:k*nij) = dest_x(1:nij)
00181 dest_y((k-1)*nij+1:k*nij) = dest_y(1:nij)
00182 dest_z((k-1)*nij+1:k*nij) = array_z(k)
00183 end do
00184
00185 else
00186
00187
00188
00189
00190 dest_x = -1
00191 dest_y = -1
00192 dest_z = -1
00193 list_entries => send_info%list_entries
00194
00195 nij = send_info%n_list
00196 nij_loc = nij - send_info%num2recv
00197
00198 #ifdef PRISM_ASSERTION
00199
00200
00201 do n = 1, nij_loc
00202 if (list_entries(1, n) < shape (1,1) .or. &
00203 list_entries(1, n) > shape (2,1)) exit
00204 end do
00205
00206 if (n < nij - send_info%num2recv) then
00207 print *, 'n, list_entry, shape', n, list_entries(1, n), &
00208 list_entries(2, n), shape
00209 call psmile_assert (__FILE__, __LINE__, &
00210 "incorrect index in list_entries")
00211 endif
00212 #endif
00213
00214 do n = 1, nij_loc
00215 dest_x (n) = min(array_x(list_entries(1,n),1),array_x(list_entries(1,n),2))
00216 dest_y (n) = min(array_y(list_entries(1,n),1),array_y(list_entries(1,n),2))
00217 dest_z (n) = array_z(list_entries(3, n))
00218 end do
00219
00220 do n = 1, nij_loc
00221 dest_x (nij+n) = max(array_x(list_entries(1,n),1),array_x(list_entries(1,n),2))
00222 dest_y (nij+n) = max(array_y(list_entries(1,n),1),array_y(list_entries(1,n),2))
00223 dest_z (nij+n) = array_z(list_entries(3, n))
00224 enddo
00225
00226 endif
00227
00228
00229
00230
00231
00232 do i = 1, nbr_cells_loc
00233 neighcells (i,2) = neighcells (i,1) + nij
00234 enddo
00235
00236 source_cell_index(:) = neighcells (:,1)
00237
00238
00239
00240
00241
00242
00243 if ( send_info%num2recv > 0 ) then
00244
00245 do n = 2, nb_corners
00246 do i = 1, send_info%num2recv
00247 neighcells (nbr_cells_loc+i,n) = &
00248 neighcells (nbr_cells_loc+i,1) + (n-1) * nij
00249 end do
00250 enddo
00251
00252
00253
00254 dest_idx_off = nij_loc
00255
00256 do l = 1, extra_search%nrecv
00257
00258 stride = nb_corners*send_info%len_sent(l)
00259
00260 do k = 1, nb_corners
00261 do n = 1, send_info%len_sent(l)
00262 idx1 = (k-1)*send_info%len_sent(l)+n
00263 idx2 = stride + idx1
00264
00265 dest_x ((k-1)*nij+dest_idx_off+n) = extra_search%dble_bufs(l)%vector(idx1)
00266 dest_y ((k-1)*nij+dest_idx_off+n) = extra_search%dble_bufs(l)%vector(idx2)
00267
00268 dest_z ((k-1)*nij+dest_idx_off+n) = array_z(1)
00269
00270 end do
00271 enddo
00272
00273
00274 dest_idx_off = dest_idx_off + send_info%len_sent(l)
00275
00276 Deallocate (extra_search%dble_bufs(l)%vector)
00277
00278 enddo
00279
00280 endif
00281
00282
00283
00284 #ifdef VERBOSE
00285 print 9980, trim(ch_id), ierror
00286
00287 call psmile_flushstd
00288 #endif /* VERBOSE */
00289
00290
00291
00292 9990 format (1x, a, ': psmile_ccompact_gauss2_dble: send_entire_valid', l2)
00293 9980 format (1x, a, ': psmile_ccompact_gauss2_dble: eof ierror = ', i8)
00294
00295 end subroutine psmile_ccompact_gauss2_dble