XBeach
|
00001 module debugging_module 00002 #if 0 00003 ! these routines need a thorough rewrite 00004 implicit none 00005 save 00006 #ifdef USEMPI 00007 interface compare 00008 module procedure comparei2 00009 module procedure comparer2 00010 module procedure comparer3 00011 end interface compare 00012 #endif 00013 00014 interface printsum 00015 module procedure printsum0 00016 module procedure printsum1 00017 module procedure printsum2 00018 module procedure printsum3 00019 module procedure printsum4 00020 module procedure printsumi0 00021 module procedure printsumi1 00022 module procedure printsumi2 00023 module procedure printsumi3 00024 module procedure printsumi4 00025 end interface printsum 00026 00027 contains 00028 00029 #ifdef USEMPI 00030 ! 00031 ! consistency check: in general, all distributed matrices and blocks 00032 ! should have the following properties: 00033 ! 00034 ! a is distributed matrix in this process 00035 ! al is distributed matrix in left neighbour process 00036 ! ar is distributed matrix in right neighbour process 00037 ! at is distributed matrix in top neighbour process 00038 ! ab is distributed matrix in bottom neighbour process 00039 ! 00040 ! a(:,2) = al(:,ny+1) 00041 ! a(2,:) = at(nx+1,:) 00042 ! except the first and last elements of these arrays 00043 ! 00044 ! 00045 ! When all these checks are ok on every process, automatically care 00046 ! has been taken for the other neighbours 00047 00048 ! 00049 ! consistency check for mnem 00050 ! if mnem = 'ALL' then all 00051 subroutine space_consistency(s,mnem,verbose) 00052 use mnemmodule 00053 use spaceparamsdef 00054 implicit none 00055 type(spacepars) :: s 00056 character(len=*) :: mnem 00057 character(len=*), intent(in), optional :: verbose 00058 integer :: j,jmin,jmax 00059 type(arraytype) :: t 00060 00061 if(mnem .eq. 'ALL') then 00062 jmin = 1 00063 jmax = numvars 00064 else 00065 jmin = chartoindex(mnem) 00066 jmax = jmin 00067 endif 00068 do j=jmin,jmax 00069 call indextos(s,j,t) 00070 select case(t%type) 00071 case('r') 00072 select case (t%rank) 00073 case(0) 00074 !call compare(t%r0,t%name) 00075 case(1) 00076 !call compare(t%r1,t%name) 00077 case(2) 00078 call compare(t%r2,t%name,verbose) 00079 case(3) 00080 call compare(t%r3,t%name,verbose) 00081 case(4) 00082 !call compare(t%r4,t%name) 00083 end select ! rank 00084 case('i') 00085 select case (t%rank) 00086 case(0) 00087 !call compare(t%i0,t%name) 00088 case(1) 00089 !call compare(t%i1,t%name) 00090 case(2) 00091 call compare(t%i2,t%name,verbose) 00092 case(3) 00093 !call compare(t%i3,t%name) 00094 case(4) 00095 !call compare(t%i4,t%name) 00096 end select ! rank 00097 end select ! type 00098 enddo 00099 00100 end subroutine space_consistency 00101 00102 subroutine comparer2(x,s,verbose) 00103 use xmpi_module 00104 use mnemmodule 00105 implicit none 00106 real*8, dimension(:,:) :: x 00107 character(len=*) :: s 00108 character(len=*),optional :: verbose 00109 00110 real*8, parameter :: eps=1.0d-60 00111 integer :: m 00112 integer :: n 00113 real*8, dimension(:,:), allocatable :: c 00114 real*8, dimension(:,:), allocatable :: r 00115 real*8, dimension(2) :: dif,difmax 00116 character*100 :: warning 00117 integer i,j 00118 00119 select case(s) 00120 case (mnem_tideinpz) 00121 return 00122 end select 00123 m=size(x,1) 00124 n=size(x,2) 00125 allocate(c(m,2)) 00126 allocate(r(2,n)) 00127 !c = x(:,1) 00128 c = x(:,1:2) 00129 !call xmpi_shift(x,':1') 00130 call xmpi_shift(x,SHIFT_Y_R,1,2) 00131 00132 !dif(1) = sum(abs(c(2:m-1)-x(2:m-1,1))) 00133 !dif(2) = sum(abs(c(1:m )-x(1:m ,1))) 00134 dif(1) = sum(abs(c(2:m-1,:)-x(2:m-1,1:2))) 00135 dif(2) = sum(abs(c(1:m ,:)-x(1:m ,1:2))) 00136 00137 call xmpi_allreduce(dif,MPI_SUM) 00138 difmax = dif 00139 00140 if(xmaster) then 00141 warning=' ' 00142 if (sum(difmax) .gt. eps) then 00143 warning = '<===++++++++++++++++++++++' 00144 endif 00145 write (*,*) 'compare (:,1:2) '//trim(s)//': ',difmax,trim(warning) 00146 endif 00147 if (present(verbose)) then 00148 if (verbose .eq. 'verbose') then 00149 if (sum(difmax) .gt. eps) then 00150 do i=1,m 00151 do j=1,2 00152 !if(abs(c(i,j)-x(i,j)) .gt. 0.0001) then 00153 if(c(i,j) .ne. 0 .or. x(i,j) .ne. 0) then 00154 print *,xmpi_rank,i,j,c(i,j),x(i,j),abs(c(i,j)-x(i,j)) 00155 endif 00156 enddo 00157 enddo 00158 endif 00159 endif 00160 endif 00161 00162 x(:,1:2) = c 00163 00164 !c = x(:,n) 00165 c = x(:,n-1:n) 00166 !call xmpi_shift(x,':n') 00167 call xmpi_shift(x,SHIFT_Y_L,3,4) 00168 00169 !dif(1) = sum(abs(c(2:m-1)-x(2:m-1,n))) 00170 !dif(2) = sum(abs(c(1:m )-x(1:m ,n))) 00171 dif(1) = sum(abs(c(2:m-1,:)-x(2:m-1,n-1:n))) 00172 dif(2) = sum(abs(c(1:m ,:)-x(1:m ,n-1:n))) 00173 !x(:,n) = c 00174 x(:,n-1:n) = c 00175 00176 call xmpi_reduce(dif,difmax,MPI_SUM) 00177 00178 if(xmaster) then 00179 warning=' ' 00180 if (sum(difmax) .gt. eps) then 00181 warning = '<===++++++++++++++++++++++' 00182 endif 00183 write (*,*) 'compare (:,n-1:n) '//trim(s)//': ',difmax,trim(warning) 00184 endif 00185 00186 !r = x(1,:) 00187 r = x(1:2,:) 00188 !call xmpi_shift(x,'1:') 00189 call xmpi_shift(x,SHIFT_X_D,1,2) 00190 00191 !dif(1) = sum(abs(r(2:n-1)-x(1,2:n-1))) 00192 !dif(2) = sum(abs(r(1:n )-x(1,1:n ))) 00193 dif(1) = sum(abs(r(:,2:n-1)-x(1:2,2:n-1))) 00194 dif(2) = sum(abs(r(:,1:n )-x(1:2,1:n ))) 00195 x(1:2,:) = r 00196 00197 call xmpi_reduce(dif,difmax,MPI_SUM) 00198 00199 if(xmaster) then 00200 warning=' ' 00201 if (sum(difmax) .gt. eps) then 00202 warning = '<===++++++++++++++++++++++' 00203 endif 00204 write (*,*) 'compare (1:2,:) '//trim(s)//': ',difmax,trim(warning) 00205 endif 00206 00207 !r = x(m,:) 00208 r = x(m-1:m,:) 00209 !call xmpi_shift(x,'m:') 00210 call xmpi_shift(x,SHIFT_X_U,3,4) 00211 00212 !dif(1) = sum(abs(r(2:n-1)-x(m,2:n-1))) 00213 !dif(2) = sum(abs(r(1:n )-x(m,1:n ))) 00214 dif(1) = sum(abs(r(:,2:n-1)-x(m-1:m,2:n-1))) 00215 dif(2) = sum(abs(r(:,1:n )-x(m-1:m,1:n ))) 00216 x(m-1:m,:) = r 00217 00218 call xmpi_reduce(dif,difmax,MPI_SUM) 00219 00220 if(xmaster) then 00221 warning=' ' 00222 if (sum(difmax) .gt. eps) then 00223 warning = '<===++++++++++++++++++++++' 00224 endif 00225 write (*,*) 'compare (m-1:m,:) '//trim(s)//': ',difmax,trim(warning) 00226 endif 00227 end subroutine comparer2 00228 00229 subroutine comparei2(x,s,verbose) 00230 use xmpi_module 00231 implicit none 00232 integer, dimension(:,:) :: x 00233 character(len=*) :: s 00234 character(len=*),optional :: verbose 00235 00236 integer, parameter :: eps=0 00237 integer :: m 00238 integer :: n 00239 integer, dimension(:), allocatable :: c 00240 integer, dimension(:), allocatable :: r 00241 integer, dimension(2) :: dif,difmax 00242 character*100 :: warning 00243 00244 print *,'comparei2 not adapted for double border scheme' 00245 return 00246 00247 if (present(verbose)) then 00248 print *,'comparei2: verbose not implemented' 00249 endif 00250 00251 m=size(x,1) 00252 n=size(x,2) 00253 allocate(c(m)) 00254 allocate(r(n)) 00255 c = x(:,1) 00256 call xmpi_shift(x,':1') 00257 00258 dif(1) = sum(abs(c(2:m-1)-x(2:m-1,1))) 00259 dif(2) = sum(abs(c(1:m )-x(1:m ,1))) 00260 x(:,1) = c 00261 00262 call xmpi_reduce(dif,difmax,MPI_SUM) 00263 00264 if(xmaster) then 00265 warning=' ' 00266 if (sum(difmax) .gt. eps) then 00267 warning = '<===++++++++++++++++++++++' 00268 endif 00269 write (*,*) 'compare (:,1) '//trim(s)//': ',difmax,trim(warning) 00270 endif 00271 00272 c = x(:,n) 00273 call xmpi_shift(x,':n') 00274 00275 dif(1) = sum(abs(c(2:m-1)-x(2:m-1,n))) 00276 dif(2) = sum(abs(c(1:m )-x(1:m ,n))) 00277 x(:,n) = c 00278 00279 call xmpi_reduce(dif,difmax,MPI_SUM) 00280 00281 if(xmaster) then 00282 warning=' ' 00283 if (sum(difmax) .gt. eps) then 00284 warning = '<===++++++++++++++++++++++' 00285 endif 00286 write (*,*) 'compare (:,n) '//trim(s)//': ',difmax,trim(warning) 00287 endif 00288 r = x(1,:) 00289 call xmpi_shift(x,'1:') 00290 00291 dif(1) = sum(abs(r(2:n-1)-x(1,2:n-1))) 00292 dif(2) = sum(abs(r(1:n )-x(1,1:n ))) 00293 x(1,:) = r 00294 00295 call xmpi_reduce(dif,difmax,MPI_SUM) 00296 00297 if(xmaster) then 00298 warning=' ' 00299 if (sum(difmax) .gt. eps) then 00300 warning = '<===++++++++++++++++++++++' 00301 endif 00302 write (*,*) 'compare (1,:) '//trim(s)//': ',difmax,trim(warning) 00303 endif 00304 00305 r = x(m,:) 00306 call xmpi_shift(x,'m:') 00307 00308 dif(1) = sum(abs(r(2:n-1)-x(m,2:n-1))) 00309 dif(2) = sum(abs(r(1:n )-x(m,1:n ))) 00310 x(m,:) = r 00311 00312 call xmpi_reduce(dif,difmax,MPI_SUM) 00313 00314 if(xmaster) then 00315 warning=' ' 00316 if (sum(difmax) .gt. eps) then 00317 warning = '<===++++++++++++++++++++++' 00318 endif 00319 write (*,*) 'compare (m,:) '//trim(s)//': ',difmax,trim(warning) 00320 endif 00321 end subroutine comparei2 00322 00323 subroutine comparer3(x,s,verbose) 00324 use xmpi_module 00325 implicit none 00326 real*8, dimension(:,:,:) :: x 00327 character(len=*) :: s 00328 character(len=*),optional :: verbose 00329 00330 real*8, parameter :: eps=1.0d-60 00331 integer :: m,n,l 00332 real*8, dimension(:,:), allocatable :: c 00333 real*8, dimension(:,:), allocatable :: r 00334 real*8, dimension(2) :: dif,difmax 00335 character*100 :: warning 00336 00337 print *,'comparei2 not adapted for double border scheme' 00338 return 00339 if (present(verbose)) then 00340 print *,'comparer3: verbose not implemented' 00341 endif 00342 m=size(x,1) 00343 n=size(x,2) 00344 l=size(x,3) 00345 00346 allocate(c(m,l)) 00347 allocate(r(n,l)) 00348 c = x(:,1,:) 00349 call xmpi_shift(x,':1') 00350 00351 dif(1) = sum(abs(c(2:m-1,:)-x(2:m-1,1,:))) 00352 dif(2) = sum(abs(c(1:m,:) -x(1:m ,1,:))) 00353 x(:,1,:) = c 00354 00355 call xmpi_reduce(dif,difmax,MPI_SUM) 00356 00357 if(xmaster) then 00358 warning=' ' 00359 if (sum(difmax) .gt. eps) then 00360 warning = '<===++++++++++++++++++++++' 00361 endif 00362 write (*,*) 'compare (:,1) '//trim(s)//': ',difmax,trim(warning) 00363 endif 00364 00365 c = x(:,n,:) 00366 call xmpi_shift(x,':n') 00367 00368 dif(1) = sum(abs(c(2:m-1,:)-x(2:m-1,n,:))) 00369 dif(2) = sum(abs(c(1:m ,:)-x(1:m ,n,:))) 00370 x(:,n,:) = c 00371 00372 call xmpi_reduce(dif,difmax,MPI_SUM) 00373 00374 if(xmaster) then 00375 warning=' ' 00376 if (sum(difmax) .gt. eps) then 00377 warning = '<===++++++++++++++++++++++' 00378 endif 00379 write (*,*) 'compare (:,n) '//trim(s)//': ',difmax,trim(warning) 00380 endif 00381 00382 r = x(1,:,:) 00383 call xmpi_shift(x,'1:') 00384 00385 dif(1) = sum(abs(r(2:n-1,:)-x(1,2:n-1,:))) 00386 dif(2) = sum(abs(r(1:n ,:)-x(1,1:n ,:))) 00387 x(1,:,:) = r 00388 00389 call xmpi_reduce(dif,difmax,MPI_SUM) 00390 00391 if(xmaster) then 00392 warning=' ' 00393 if (sum(difmax) .gt. eps) then 00394 warning = '<===++++++++++++++++++++++' 00395 endif 00396 write (*,*) 'compare (1,:) '//trim(s)//': ',difmax,trim(warning) 00397 endif 00398 00399 r = x(m,:,:) 00400 call xmpi_shift(x,'m:') 00401 00402 dif(1) = sum(abs(r(2:n-1,:)-x(m,2:n-1,:))) 00403 dif(2) = sum(abs(r(1:n ,:)-x(m,1:n ,:))) 00404 x(m,:,:) = r 00405 00406 call xmpi_reduce(dif,difmax,MPI_SUM) 00407 00408 if(xmaster) then 00409 warning=' ' 00410 if (sum(difmax) .gt. eps) then 00411 warning = '<===++++++++++++++++++++++' 00412 endif 00413 write (*,*) 'compare (m,:) '//trim(s)//': ',difmax,trim(warning) 00414 endif 00415 end subroutine comparer3 00416 ! 00417 #endif 00418 00419 subroutine printsum0(f,str,id,val) 00420 implicit none 00421 integer, intent(in) :: f 00422 character(*),intent(in) :: str 00423 integer, intent(in) :: id 00424 real*8, intent(in) :: val 00425 write(f,*) 'printsum ',id,' ',str,':',val 00426 end subroutine printsum0 00427 00428 subroutine printsum1(f,str,id,val) 00429 implicit none 00430 integer, intent(in) :: f 00431 character(*), intent(in) :: str 00432 integer, intent(in) :: id 00433 real*8, pointer, dimension(:) :: val 00434 if (associated(val) ) then 00435 write(f,*) 'printsum ',id,' ',str,':',sum(val),shape(val) 00436 else 00437 write(f,*) 'printsum ',id,' ',str,':',' Not allocated' 00438 endif 00439 end subroutine printsum1 00440 00441 subroutine printsum2(f,str,id,val) 00442 implicit none 00443 integer, intent(in) :: f 00444 character(*), intent(in) :: str 00445 integer, intent(in) :: id 00446 real*8, pointer, dimension(:,:) :: val 00447 if (associated(val) ) then 00448 write(f,*) 'printsum ',id,' ',str,':',sum(val),shape(val) 00449 else 00450 write(f,*) 'printsum ',id,' ',str,':',' Not allocated' 00451 endif 00452 00453 end subroutine printsum2 00454 00455 subroutine printsum3(f,str,id,val) 00456 implicit none 00457 integer, intent(in) :: f 00458 character(*), intent(in) :: str 00459 integer, intent(in) :: id 00460 real*8, pointer, dimension(:,:,:) :: val 00461 if (associated(val) ) then 00462 write(f,*) 'printsum ',id,' ',str,':',sum(val),shape(val) 00463 else 00464 write(f,*) 'printsum ',id,' ',str,':',' Not allocated' 00465 endif 00466 end subroutine printsum3 00467 00468 subroutine printsum4(f,str,id,val) 00469 implicit none 00470 integer, intent(in) :: f 00471 character(*), intent(in) :: str 00472 integer, intent(in) :: id 00473 real*8, pointer, dimension(:,:,:,:) :: val 00474 if (associated(val) ) then 00475 write(f,*) 'printsum ',id,' ',str,':',sum(val),shape(val) 00476 else 00477 write(f,*) 'printsum ',id,' ',str,':',' Not allocated' 00478 endif 00479 end subroutine printsum4 00480 00481 subroutine printsumi0(f,str,id,val) 00482 implicit none 00483 integer, intent(in) :: f 00484 character(*),intent(in) :: str 00485 integer, intent(in) :: id 00486 integer, intent(in) :: val 00487 write(f,*) 'printsum ',id,' ',str,':',val 00488 end subroutine printsumi0 00489 00490 subroutine printsumi1(f,str,id,val) 00491 implicit none 00492 integer, intent(in) :: f 00493 character(*), intent(in) :: str 00494 integer, intent(in) :: id 00495 integer, pointer, dimension(:) :: val 00496 if (associated(val) ) then 00497 write(f,*) 'printsum ',id,' ',str,':',sum(val),shape(val) 00498 else 00499 write(f,*) 'printsum ',id,' ',str,':',' Not allocated' 00500 endif 00501 end subroutine printsumi1 00502 00503 subroutine printsumi2(f,str,id,val) 00504 implicit none 00505 integer, intent(in) :: f 00506 character(*), intent(in) :: str 00507 integer, intent(in) :: id 00508 integer, pointer, dimension(:,:) :: val 00509 if (associated(val) ) then 00510 write(f,*) 'printsum ',id,' ',str,':',sum(val),shape(val) 00511 else 00512 write(f,*) 'printsum ',id,' ',str,':',' Not allocated' 00513 endif 00514 00515 end subroutine printsumi2 00516 00517 subroutine printsumi3(f,str,id,val) 00518 implicit none 00519 integer, intent(in) :: f 00520 character(*), intent(in) :: str 00521 integer, intent(in) :: id 00522 integer, pointer, dimension(:,:,:) :: val 00523 if (associated(val) ) then 00524 write(f,*) 'printsum ',id,' ',str,':',sum(val),shape(val) 00525 else 00526 write(f,*) 'printsum ',id,' ',str,':',' Not allocated' 00527 endif 00528 end subroutine printsumi3 00529 00530 subroutine printsumi4(f,str,id,val) 00531 implicit none 00532 integer, intent(in) :: f 00533 character(*), intent(in) :: str 00534 integer, intent(in) :: id 00535 integer, allocatable, dimension(:,:,:,:) :: val 00536 if (allocated(val) ) then 00537 write(f,*) 'printsum ',id,' ',str,':',sum(val),shape(val) 00538 else 00539 write(f,*) 'printsum ',id,' ',str,':',' Not allocated' 00540 endif 00541 end subroutine printsumi4 00542 00543 subroutine printssums(s,str) 00544 #ifdef USEMPI 00545 use xmpi_module 00546 #endif 00547 use mnemmodule 00548 use spaceparams 00549 type (spacepars), intent(in) :: s 00550 character(*), intent(in) :: str 00551 integer :: id, f,i 00552 type(arraytype) :: t 00553 00554 #ifdef USEMPI 00555 id = xmpi_rank 00556 if (id .gt. 0 ) then 00557 return 00558 endif 00559 #else 00560 id=0 00561 #endif 00562 00563 f = 100+4*numvars+id+1 00564 00565 write(f, *) 'printsum: ',id,'Start of printssums ',str 00566 do i=1,numvars 00567 call indextos(s,i,t) 00568 select case (t%rank) 00569 case (0) 00570 if (t%type .eq. 'i') then 00571 call printsum(f,t%name,id,t%i0) 00572 else 00573 call printsum(f,t%name,id,t%r0) 00574 endif 00575 case (1) 00576 if (t%type .eq. 'i') then 00577 call printsum(f,t%name,id,t%i1) 00578 else 00579 call printsum(f,t%name,id,t%r1) 00580 endif 00581 case (2) 00582 if (t%type .eq. 'i') then 00583 call printsum(f,t%name,id,t%i2) 00584 else 00585 call printsum(f,t%name,id,t%r2) 00586 endif 00587 case (3) 00588 if (t%type .eq. 'i') then 00589 call printsum(f,t%name,id,t%i3) 00590 else 00591 call printsum(f,t%name,id,t%r3) 00592 endif 00593 case (4) 00594 if (t%type .eq. 'i') then 00595 call printsum(f,t%name,id,t%i4) 00596 else 00597 call printsum(f,t%name,id,t%r4) 00598 endif 00599 end select 00600 enddo 00601 #ifdef USEMPI 00602 call printsum(f,'s%is',id,s%is) 00603 call printsum(f,'s%js',id,s%js) 00604 call printsum(f,'s%lm',id,s%lm) 00605 call printsum(f,'s%ln',id,s%ln) 00606 #endif 00607 end subroutine printssums 00608 00609 subroutine printssumso(s) 00610 use spaceparams 00611 implicit none 00612 type (spacepars), intent(in) :: s 00613 00614 write(*,*)'Start of printssumso' 00615 00616 write(*,*)'s%xz',sum(s%xz) 00617 write(*,*)'s%yz',sum(s%yz) 00618 write(*,*)'s%zs',sum(s%zs) 00619 write(*,*)'s%u',sum(s%u) 00620 write(*,*)'s%v',sum(s%v) 00621 write(*,*)'s%ue',sum(s%ue) 00622 write(*,*)'s%ve',sum(s%ve) 00623 write(*,*)'s%H',sum(s%H) 00624 write(*,*)'s%urms',sum(s%urms) 00625 write(*,*)'s%zb',sum(s%zb) 00626 write(*,*)'s%hh',sum(s%hh) 00627 write(*,*)'s%Fx',sum(s%Fx) 00628 write(*,*)'s%Fy',sum(s%Fy) 00629 write(*,*)'s%E',sum(s%E) 00630 write(*,*)'s%R',sum(s%R) 00631 write(*,*)'s%D',sum(s%D) 00632 end subroutine printssumso 00633 #endif 00634 end module debugging_module