1! --- 2! Copyright (C) 1996-2016 The SIESTA group 3! This file is distributed under the terms of the 4! GNU General Public License: see COPYING in the top directory 5! or http://www.gnu.org/copyleft/gpl.txt . 6! See Docs/Contributors.txt for a list of contributors. 7! --- 8!!@LICENSE 9! 10! ================================================================== 11! Allocation, reallocation, and deallocation utility routines 12! Written by J.M.Soler. May 2000. 13! ================================================================== 14! SUBROUTINE alloc_default( old, new, restore, & 15! copy, shrink, imin, routine ) 16! Sets defaults for allocation 17! INPUT (optional): 18! type(allocDefaults) restore : default settings to be restored 19! logical copy : Copy old array to new array? 20! logical shrink : Reduce array size? 21! integer imin : First index (typically 1 in Fortan, 22! 0 in C) 23! character(len=*) routine : Name of calling routine 24! OUTPUT (optional): 25! type(allocDefaults) old : default settings before the call 26! type(allocDefaults) new : default settings after the call 27! BEHAVIOR: 28! All these defaults can be superseeded by optional arguments in 29! each call to re_alloc. 30! Initial default values: copy = .true. 31! shrink = .true. 32! imin = 1 33! routine = 'unknown' 34! If restore is present together with any of copy, shrink, imin, or 35! routine, these are applied AFTER resetting the restore defaults. 36! USAGE: 37! In order to restore the allocation defaults possibly set by the 38! calling routine, the suggested construction is: 39! use alloc_module 40! type(allocDefaults) oldDefaults 41! call alloc_default( old=oldDefaults, routine=..., & 42! copy=..., shrink=... ) 43! call re_alloc(...) 44! call alloc_default( restore=oldDefaults ) 45! Notice that, if the restore call is skipped, the new defaults will 46! stay in effect until a new call to alloc_dafault is made. 47! ================================================================== 48! SUBROUTINE alloc_report( level, unit, file, printNow, threshold, shutdown ) 49! Sets the output file for the allocation report 50! INPUT (optional): 51! integer :: level : Level (detail) of report 52! integer :: unit : Output file unit 53! character*(*):: file : Output file name 54! logical :: printNow : If present & true => print report now 55! real(dp) :: threshold : Memory threshold (in bytes) to print 56! the memory use of any given array 57! logical :: shutdown : If present & true then close output unit. 58! BEHAVIOR: 59! The detail/extent of the report increses with the value of level: 60! level=0 : no report at all (the default) 61! level=1 : only total memory peak and where it occurred 62! level=2 : detailed report created but printed only upon request 63! level=3 : detailed report printed at every new memory peak 64! level=4 : print every individual reallocation or deallocation 65! If unit is present, alloc_report merely takes note of it for 66! future use, assuming that it has been already open outside. 67! In this case, file is not used. 68! If unit is absent, and file is present, a file with that 69! name is open for future use. 70! If both arguments are absent, a file named 'alloc_report' 71! is open for future use. 72! If alloc_report is called with printNow=.true. several times in 73! a program, with the same unit or file argument, the subsequent 74! reports are written consecutively in the same file, each with a 75! time stamp header. 76! If threshold is not present, threshold=0 is assumed. 77! When performing library invocations of this module in say 78! external methods it is vital to close units before exiting the 79! library to ensure that no dangling files are kept open. 80! To ensure the file is closed before return one can use 81! the shutdown=.true. argument to forcefully close (and free) 82! the unit used to report. If not present or .false. the unit 83! will be kept open. 84! In parallel execution, the report sections that involve every 85! reallocation (levels 1, 3, and 4) are written only by node 0. 86! The section that is written upon request (level 2) is written 87! only by the node with the highest peak of memory up to that time, 88! but it contains a summary of the memory used by all other nodes. 89! In parallel execution, the nodes that share the same file 90! system (e.g. different chip cores or NFS-connected nodes) write 91! on the same file. Otherwise they write on files with the same name 92! in their local disks. 93! ================================================================== 94! SUBROUTINE re_alloc( array, [i1min,] i1max, 95! [[i2min,] i2max, [[i3min,] i3max]], 96! name, routine, copy, shrink ) 97! INPUT: 98! integer :: i1min : Lower bound of first dimension 99! If not present, it is fixed by 100! the last call to alloc_default. 101! If present and the rank is 2(3), 102! then i2min(&i3min) must also be 103! present 104! integer :: i1max : Upper bound of first dimension 105! integer :: i2min,i2max : Bounds of second dimension, if 106! applicable 107! integer :: i3min,i3max : Bounds of third dimension, if appl. 108! 109! INPUT (optional): 110! character*(*) :: name : Actual array name or a label for it 111! character*(*) :: routine : Name of the calling routine 112! or routine section 113! logical :: copy : Save (copy) contents of old array 114! to new array? 115! logical :: shrink : Reallocate if the new array bounds 116! are contained within the old ones? 117! If not present, copy and/or shrink 118! are fixed by the last call to 119! alloc_default. 120! INPUT/OUTPUT: 121! TYPE, pointer :: array : Array to be allocated or reallocated. 122! Implemented types and ranks are: 123! integer, rank 1, 2, 3 124! integer*8, rank 1 125! real*4, rank 1, 2, 3, 4 126! real*8, rank 1, 2, 3, 4 127! complex*16, rank 1, 2 128! logical, rank 1, 2, 3 129! character(len=*), rank 1 130! BEHAVIOR: 131! Pointers MUST NOT enter in an undefined state. Before using them 132! for the first time, they must be nullified explicitly. Alternatively, 133! in f95, they can be initialized as null() upon declaration. 134! If argument array is not associated on input, it is just allocated. 135! If array is associated and has the same bounds (or smaller bonds 136! and shrink is false) nothing is done. Thus, it is perfectly safe and 137! efficient to call re_alloc repeatedly without deallocating the array. 138! However, subroutine dealloc is provided to eliminate large arrays 139! when they are not needed. 140! In order to save (copy) the contents of the old array, the new array 141! needs to be allocated before deallocating the old one. Thus, if the 142! contents are not needed, or if reducing memory is a must, calling 143! re_alloc with copy=.false. makes it to deallocate before allocating. 144! The elements that are not copied (because copy=.false. or because 145! they are outside the bounds of the input array) return with value 146! zero (integer and real), .false. (logical), or blank (character). 147! If imin>imax for any dimension, the array pointer returns 148! associated to a zero-size array. 149! Besides allocating or reallocating the array, re_alloc keeps a count 150! of memory usage (and prints a report in a file previously fixed by a 151! call to alloc_report). Thus, it is not recommended to call re_alloc 152! to reallocate an array not allocated by it the first time. 153! If name is not present, the memory count associated to the 154! allocation/deallocation is still made, but the allocation report 155! will account the array under the generic name 'unknown'. 156! The routine argument is NOT used to classify the counting of 157! memory usage. The classification uses only the name argument. 158! This is because a pointer may be allocated in one routine and 159! deallocated in a different one (i.e. when it is used to return an 160! array whose size is not known in advance). However, the routine 161! argument is reported when alloc_report=4, and it is also used to 162! report in which routine the memory peack occurs. If you want the 163! routine name to be used for classification, you should include it 164! as part of the name argument, like in name='matInvert '//'aux'. 165! ==================================================================--- 166! SUBROUTINE de_alloc( array, name, routine ) 167! INPUT (optional): 168! character*(*) :: name : Actual array name or a label for it 169! character*(*) :: routine : Name of the calling routine 170! or routine section 171! INPUT/OUTPUT: 172! TYPE, pointer :: array : Array be deallocated (same types and 173! kinds as in re_alloc). 174! BEHAVIOR: 175! Besides deallocating the array, re_alloc decreases the count of 176! memory usage previously counted by re_alloc. Thus, dealloc should 177! not be called to deallocate an array not allocated by re_alloc. 178! Equally, arrays allocated or reallocated by re_alloc should be 179! deallocated by dealloc. 180! ==================================================================--- 181! SUBROUTINE alloc_count( delta_size, type, name, routine ) 182! INPUT: 183! integer :: delta_size : +/-size(array) 184! + => allocation 185! - => deallocation 186! character :: type : 'I' => integer 187! 'E' => integer*8 188! 'R' => real*4 189! 'D' => real*8 190! 'L' => logical 191! 'S' => character (string) 192! INPUT (optional): 193! character(len=*) :: name : Actual array name or a label for it 194! character(len=*) :: routine : Name of the calling routine 195! or routine section 196! USAGE: 197! integer, allocatable:: intArray(:) 198! double precision,allocatable:: doubleArray(:) 199! complex, allocatable:: complexArray(:) 200! allocate( intArray(n), doubleArray(n), complexArray(n) ) 201! call alloc_count( +n, 'I', 'intArray', programName ) 202! call alloc_count( +n, 'D', 'doubleArray', programName ) 203! call alloc_count( +2*n, 'R', 'complexArray', programName ) 204! deallocate( intArray, doubleArray, complexArray ) 205! call alloc_count( -n, 'I', 'intArray', programName ) 206! call alloc_count( -n, 'D', 'doubleArray', programName ) 207! call alloc_count( -2*n, 'R', 'complexArray', programName ) 208! 209! ==================================================================--- 210 211MODULE alloc 212 213 use precision, only: sp ! Single precision real type 214 use precision, only: dp ! Double precision real type 215 use precision, only: i8b ! 8-byte integer 216 use parallel, only: Node ! My processor node index 217 use parallel, only: Nodes ! Number of parallel processors 218 use parallel, only: ionode ! Am I the I/O processor? 219 use parallel, only: parallel_init ! Initialize parallel variables 220 use sys, only: die ! Termination routine 221 use m_io, only: io_assign ! Get and reserve an available IO unit 222#ifdef MPI 223 use mpi_siesta 224#endif 225 226 implicit none 227 228PUBLIC :: & 229 alloc_default, &! Sets allocation defaults 230 alloc_report, &! Sets log report defaults 231 re_alloc, &! Allocation/reallocation 232 de_alloc, &! Deallocation 233 alloc_count, &! Memory counting for external allocs 234 allocDefaults ! Derived type to hold allocation defaults 235 236PRIVATE ! Nothing is declared public beyond this point 237 238 interface de_alloc 239 module procedure & 240 dealloc_i1, dealloc_i2, dealloc_i3, & 241 dealloc_E1, & 242 dealloc_r1, dealloc_r2, dealloc_r3, dealloc_r4, & 243 dealloc_d1, dealloc_d2, dealloc_d3, dealloc_d4, & 244 dealloc_c1, dealloc_c2, & 245 dealloc_z1, dealloc_z2, dealloc_z3, dealloc_z4, & 246 dealloc_l1, dealloc_l2, dealloc_l3, & 247 dealloc_s1 248 end interface 249 250 interface re_alloc 251 module procedure & 252 realloc_i1, realloc_i2, realloc_i3, & 253 realloc_E1, & 254 realloc_r1, realloc_r2, realloc_r3, realloc_r4, & 255 realloc_d1, realloc_d2, realloc_d3, realloc_d4, & 256 realloc_c1, realloc_c2, & 257 realloc_z1, realloc_z2, realloc_z3, realloc_z4, & 258 realloc_l1, realloc_l2, realloc_l3, & 259 realloc_s1 260! module procedure & ! AG: Dangerous!!! 261! realloc_i1s, realloc_i2s, realloc_i3s, & 262! realloc_r1s, realloc_r2s, realloc_r3s, realloc_r4s, & 263! realloc_d1s, realloc_d2s, realloc_d3s, realloc_d4s, & 264! realloc_l1s, realloc_l2s, realloc_l3s 265 end interface 266 267 ! Initial default values 268 character(len=*), parameter :: & 269 DEFAULT_NAME = 'unknown_name' ! Array name default 270 character(len=*), parameter :: & 271 DEFAULT_ROUTINE = 'unknown_routine' ! Routine name default 272 integer, save :: & 273 REPORT_LEVEL = 0, &! Level (detail) of allocation report 274 REPORT_UNIT = 0 ! Output file unit for report 275 character(len=50), save :: & 276 REPORT_FILE = 'alloc_report' ! Output file name for report 277 real(dp), save :: & 278 REPORT_THRESHOLD = 0 ! Memory threshold (in bytes) to print 279 ! the memory use of any given array 280 281 ! Derived type to hold allocation default options 282 type allocDefaults 283 private 284 logical :: copy = .true. ! Copy default 285 logical :: shrink = .true. ! Shrink default 286 integer :: imin = 1 ! Imin default 287 character(len=32):: routine = DEFAULT_ROUTINE ! Routine name default 288 end type allocDefaults 289 290 ! Object to hold present allocation default options 291 type(allocDefaults), save :: DEFAULT 292 293 ! Internal auxiliary type for a binary tree 294 type TREE 295 character(len=80) :: name ! Name of an allocated array 296 real(DP) :: mem ! Present memory use of the array 297 real(DP) :: max ! Maximum memory use of the array 298 real(DP) :: peak ! Memory use of the array during 299 ! peak of total memory 300 type(TREE), pointer :: left ! Pointer to data of allocated arrays 301 ! preceeding in alphabetical order 302 type(TREE), pointer :: right ! Pointer to data of allocated arrays 303 ! trailing in alphabetical order 304 end type TREE 305 306 ! Global variables used to store allocation data 307 real(DP), parameter :: MBYTE = 1.e6_dp 308 type(TREE), pointer, save :: REPORT_TREE 309 real(DP), save :: TOT_MEM = 0._dp 310 real(DP), save :: PEAK_MEM = 0._dp 311 character(len=80), save :: PEAK_ARRAY = ' ' 312 character(len=32), save :: PEAK_ROUTINE = ' ' 313 integer, save :: MAX_LEN = 0 314 315 ! Other common variables 316 integer :: IERR 317 logical :: ASSOCIATED_ARRAY, NEEDS_ALLOC, NEEDS_COPY, NEEDS_DEALLOC 318 319CONTAINS 320 321! ================================================================== 322 323SUBROUTINE alloc_default( old, new, restore, & 324 routine, copy, shrink, imin ) 325implicit none 326type(allocDefaults), optional, intent(out) :: old, new 327type(allocDefaults), optional, intent(in) :: restore 328character(len=*), optional, intent(in) :: routine 329logical, optional, intent(in) :: copy, shrink 330integer, optional, intent(in) :: imin 331 332if (present(old)) old = DEFAULT 333if (present(restore)) DEFAULT = restore 334if (present(copy)) DEFAULT%copy = copy 335if (present(shrink)) DEFAULT%shrink = shrink 336if (present(imin)) DEFAULT%imin = imin 337if (present(routine)) DEFAULT%routine = routine 338if (present(new)) new = DEFAULT 339 340END SUBROUTINE alloc_default 341 342! ================================================================== 343 344SUBROUTINE alloc_report( level, unit, file, printNow, threshold, & 345 shutdown) 346 347implicit none 348 349integer, optional, intent(in) :: level, unit 350character(len=*), optional, intent(in) :: file 351logical, optional, intent(in) :: printNow 352real(dp), optional, intent(in) :: threshold 353logical, optional, intent(in) :: shutdown 354 355logical :: is_open 356 357#ifdef MPI 358integer :: MPIerror 359#endif 360 361if (present(level)) then 362 REPORT_LEVEL = level 363end if 364 365if (REPORT_LEVEL <= 0) RETURN 366 367if (node == 0) then 368 if (present(unit)) then ! Assume that unit has been open outside 369 if (unit > 0) then 370 REPORT_UNIT = unit 371 if (present(file)) then 372 REPORT_FILE = file 373 else 374 REPORT_FILE = 'unknown' 375 end if 376 end if 377 else if (present(file)) then ! If file is the same, do nothing 378 if (file /= REPORT_FILE) then ! Check if file was open outside 379 REPORT_FILE = file 380 inquire( file=REPORT_FILE, opened=is_open, number=REPORT_UNIT ) 381 if (.not.is_open) then ! Open new file 382 call io_assign(REPORT_UNIT) 383 open( REPORT_UNIT, file=REPORT_FILE, status='unknown') 384 write(REPORT_UNIT,*) ' ' ! Overwrite previous reports 385 end if 386 end if 387 else if (REPORT_UNIT==0) then ! No unit has been open yet 388 REPORT_FILE = 'alloc_report' 389 call io_assign(REPORT_UNIT) 390 open( REPORT_UNIT, file=REPORT_FILE, status='unknown') 391 write(REPORT_UNIT,*) ' ' ! Overwrite previous reports 392 end if 393end if 394 395#ifdef MPI 396! Distribute information to other nodes and open REPORT_UNIT 397! NP: 398! I am not too happy about this 399! This forces a certain unit to be Bcasted to the others. 400! If that unit is already open on the other nodes then you will have 401! this module and some other module write to the same file. 402! Perhaps we should just open the file from this node with 403! the node id appended? 404call MPI_Bcast(REPORT_FILE,50,MPI_character,0,MPI_Comm_World,MPIerror) 405! JMS: open file only in node 0 406!if (node > 0) then 407! open( REPORT_UNIT, file=REPORT_FILE ) 408!end if 409#endif 410 411if (present(threshold)) REPORT_THRESHOLD = threshold 412 413if (present(printNow)) then 414 if (printNow) call print_report( ) 415end if 416 417if (present(shutdown)) then 418 if ( shutdown .and. REPORT_UNIT /= 0 ) then 419 inquire( REPORT_UNIT, opened=is_open ) 420 if ( is_open ) call io_close(REPORT_UNIT) 421 end if 422end if 423 424END SUBROUTINE alloc_report 425 426! ================================================================== 427! Integer array reallocs 428! ================================================================== 429 430SUBROUTINE realloc_i1( array, i1min, i1max, & 431 name, routine, copy, shrink ) 432! Arguments 433implicit none 434integer, dimension(:), pointer :: array 435integer, intent(in) :: i1min 436integer, intent(in) :: i1max 437character(len=*), optional, intent(in) :: name 438character(len=*), optional, intent(in) :: routine 439logical, optional, intent(in) :: copy 440logical, optional, intent(in) :: shrink 441 442! Internal variables and arrays 443character, parameter :: type='I' 444integer, parameter :: rank=1 445integer, dimension(:), pointer :: old_array 446integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 447 448! Get old array bounds 449ASSOCIATED_ARRAY = associated(array) 450if (ASSOCIATED_ARRAY) then 451 old_array => array ! Keep pointer to old array 452 old_bounds(1,:) = lbound(old_array) 453 old_bounds(2,:) = ubound(old_array) 454end if 455 456! Copy new requested array bounds 457new_bounds(1,:) = (/ i1min /) 458new_bounds(2,:) = (/ i1max /) 459 460! Find if it is a new allocation or a true reallocation, 461! and if the contents need to be copied (saved) 462! Argument b returns common bounds 463! Options routine also reads common variable ASSOCIATED_ARRAY, 464! and it sets NEEDS_ALLOC, NEEDS_DEALLOC, and NEEDS_COPY 465call options( b, c, old_bounds, new_bounds, copy, shrink ) 466! Deallocate old space 467if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 468 call alloc_count( -size(old_array), type, name, routine ) 469 deallocate(old_array,stat=IERR) 470 call alloc_err( IERR, name, routine, old_bounds ) 471end if 472 473! Allocate new space 474if (NEEDS_ALLOC) then 475 allocate( array(b(1,1):b(2,1)), stat=IERR ) 476 call alloc_err( IERR, name, routine, new_bounds ) 477 call alloc_count( size(array), type, name, routine ) 478 array = 0 479end if 480 481! Copy contents and deallocate old space 482if (NEEDS_COPY) then 483 array(c(1,1):c(2,1)) = old_array(c(1,1):c(2,1)) 484 call alloc_count( -size(old_array), type, name, routine ) 485 deallocate(old_array,stat=IERR) 486 call alloc_err( IERR, name, routine, old_bounds ) 487end if 488END SUBROUTINE realloc_i1 489 490! ================================================================== 491SUBROUTINE realloc_i2( array, i1min,i1max, i2min,i2max, & 492 name, routine, copy, shrink ) 493implicit none 494character, parameter :: type='I' 495integer, parameter :: rank=2 496integer, dimension(:,:), pointer :: array, old_array 497integer, intent(in) :: i1min, i1max, i2min, i2max 498character(len=*), optional, intent(in) :: name, routine 499logical, optional, intent(in) :: copy, shrink 500integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 501ASSOCIATED_ARRAY = associated(array) 502if (ASSOCIATED_ARRAY) then 503 old_array => array 504 old_bounds(1,:) = lbound(old_array) 505 old_bounds(2,:) = ubound(old_array) 506end if 507new_bounds(1,:) = (/ i1min, i2min /) 508new_bounds(2,:) = (/ i1max, i2max /) 509call options( b, c, old_bounds, new_bounds, copy, shrink ) 510if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 511 call alloc_count( -size(old_array), type, name, routine ) 512 deallocate(old_array,stat=IERR) 513 call alloc_err( IERR, name, routine, old_bounds ) 514end if 515if (NEEDS_ALLOC) then 516 allocate( array(b(1,1):b(2,1),b(1,2):b(2,2)), stat=IERR ) 517 call alloc_err( IERR, name, routine, new_bounds ) 518 call alloc_count( size(array), type, name, routine ) 519 array = 0 520end if 521if (NEEDS_COPY) then 522 array(c(1,1):c(2,1),c(1,2):c(2,2)) = & 523 old_array(c(1,1):c(2,1),c(1,2):c(2,2)) 524 call alloc_count( -size(old_array), type, name, routine ) 525 deallocate(old_array,stat=IERR) 526 call alloc_err( IERR, name, routine, old_bounds ) 527end if 528END SUBROUTINE realloc_i2 529! ================================================================== 530 531SUBROUTINE realloc_i3( array, i1min,i1max, i2min,i2max, i3min,i3max, & 532 name, routine, copy, shrink ) 533implicit none 534character, parameter :: type='I' 535integer, parameter :: rank=3 536integer, dimension(:,:,:), pointer :: array, old_array 537integer, intent(in) :: i1min,i1max, i2min,i2max, & 538 i3min,i3max 539character(len=*), optional, intent(in) :: name, routine 540logical, optional, intent(in) :: copy, shrink 541integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 542ASSOCIATED_ARRAY = associated(array) 543if (ASSOCIATED_ARRAY) then 544 old_array => array 545 old_bounds(1,:) = lbound(old_array) 546 old_bounds(2,:) = ubound(old_array) 547end if 548new_bounds(1,:) = (/ i1min, i2min, i3min /) 549new_bounds(2,:) = (/ i1max, i2max, i3max /) 550call options( b, c, old_bounds, new_bounds, copy, shrink ) 551if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 552 call alloc_count( -size(old_array), type, name, routine ) 553 deallocate(old_array,stat=IERR) 554 call alloc_err( IERR, name, routine, old_bounds ) 555end if 556if (NEEDS_ALLOC) then 557 allocate( array(b(1,1):b(2,1),b(1,2):b(2,2),b(1,3):b(2,3)),stat=IERR) 558 call alloc_err( IERR, name, routine, new_bounds ) 559 call alloc_count( size(array), type, name, routine ) 560 array = 0 561end if 562if (NEEDS_COPY) then 563 array(c(1,1):c(2,1),c(1,2):c(2,2),c(1,3):c(2,3)) = & 564 old_array(c(1,1):c(2,1),c(1,2):c(2,2),c(1,3):c(2,3)) 565 call alloc_count( -size(old_array), type, name, routine ) 566 deallocate(old_array,stat=IERR) 567 call alloc_err( IERR, name, routine, old_bounds ) 568end if 569END SUBROUTINE realloc_i3 570! ================================================================== 571SUBROUTINE realloc_E1( array, i1min, i1max, & 572 name, routine, copy, shrink ) 573! Arguments 574implicit none 575integer(i8b), dimension(:), pointer :: array 576integer, intent(in) :: i1min 577integer, intent(in) :: i1max 578character(len=*), optional, intent(in) :: name 579character(len=*), optional, intent(in) :: routine 580logical, optional, intent(in) :: copy 581logical, optional, intent(in) :: shrink 582 583! Internal variables and arrays 584character, parameter :: type='I' 585integer, parameter :: rank=1 586integer(i8b), dimension(:), pointer :: old_array 587integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 588 589! Get old array bounds 590ASSOCIATED_ARRAY = associated(array) 591if (ASSOCIATED_ARRAY) then 592 old_array => array ! Keep pointer to old array 593 old_bounds(1,:) = lbound(old_array) 594 old_bounds(2,:) = ubound(old_array) 595end if 596 597! Copy new requested array bounds 598new_bounds(1,:) = (/ i1min /) 599new_bounds(2,:) = (/ i1max /) 600 601! Find if it is a new allocation or a true reallocation, 602! and if the contents need to be copied (saved) 603! Argument b returns common bounds 604! Options routine also reads common variable ASSOCIATED_ARRAY, 605! and it sets NEEDS_ALLOC, NEEDS_DEALLOC, and NEEDS_COPY 606call options( b, c, old_bounds, new_bounds, copy, shrink ) 607 608! Deallocate old space 609if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 610 call alloc_count( -size(old_array), type, name, routine ) 611 deallocate(old_array,stat=IERR) 612 call alloc_err( IERR, name, routine, old_bounds ) 613end if 614 615! Allocate new space 616if (NEEDS_ALLOC) then 617 allocate( array(b(1,1):b(2,1)), stat=IERR ) 618 call alloc_err( IERR, name, routine, new_bounds ) 619 call alloc_count( size(array), type, name, routine ) 620 array = 0 621end if 622 623! Copy contents and deallocate old space 624if (NEEDS_COPY) then 625 array(c(1,1):c(2,1)) = old_array(c(1,1):c(2,1)) 626 call alloc_count( -size(old_array), type, name, routine ) 627 deallocate(old_array,stat=IERR) 628 call alloc_err( IERR, name, routine, old_bounds ) 629end if 630 631END SUBROUTINE realloc_E1 632 633! ================================================================== 634! Single precision real array reallocs 635! ================================================================== 636SUBROUTINE realloc_r1( array, i1min, i1max, & 637 name, routine, copy, shrink ) 638implicit none 639character, parameter :: type='R' 640integer, parameter :: rank=1 641real(SP), dimension(:), pointer :: array, old_array 642integer, intent(in) :: i1min, i1max 643character(len=*), optional, intent(in) :: name, routine 644logical, optional, intent(in) :: copy, shrink 645integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 646ASSOCIATED_ARRAY = associated(array) 647if (ASSOCIATED_ARRAY) then 648 old_array => array 649 old_bounds(1,:) = lbound(old_array) 650 old_bounds(2,:) = ubound(old_array) 651end if 652new_bounds(1,:) = (/ i1min /) 653new_bounds(2,:) = (/ i1max /) 654call options( b, c, old_bounds, new_bounds, copy, shrink ) 655if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 656 call alloc_count( -size(old_array), type, name, routine ) 657 deallocate(old_array,stat=IERR) 658 call alloc_err( IERR, name, routine, old_bounds ) 659end if 660if (NEEDS_ALLOC) then 661 allocate( array(b(1,1):b(2,1)), stat=IERR ) 662 call alloc_err( IERR, name, routine, new_bounds ) 663 call alloc_count( size(array), type, name, routine ) 664 array = 0._sp 665end if 666if (NEEDS_COPY) then 667 array(c(1,1):c(2,1)) = old_array(c(1,1):c(2,1)) 668 call alloc_count( -size(old_array), type, name, routine ) 669 deallocate(old_array,stat=IERR) 670 call alloc_err( IERR, name, routine, old_bounds ) 671end if 672END SUBROUTINE realloc_r1 673! ================================================================== 674SUBROUTINE realloc_r2( array, i1min,i1max, i2min,i2max, & 675 name, routine, copy, shrink ) 676implicit none 677character, parameter :: type='R' 678integer, parameter :: rank=2 679real(SP), dimension(:,:), pointer :: array, old_array 680integer, intent(in) :: i1min, i1max, i2min, i2max 681character(len=*), optional, intent(in) :: name, routine 682logical, optional, intent(in) :: copy, shrink 683integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 684ASSOCIATED_ARRAY = associated(array) 685if (ASSOCIATED_ARRAY) then 686 old_array => array 687 old_bounds(1,:) = lbound(old_array) 688 old_bounds(2,:) = ubound(old_array) 689end if 690new_bounds(1,:) = (/ i1min, i2min /) 691new_bounds(2,:) = (/ i1max, i2max /) 692call options( b, c, old_bounds, new_bounds, copy, shrink ) 693if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 694 call alloc_count( -size(old_array), type, name, routine ) 695 deallocate(old_array,stat=IERR) 696 call alloc_err( IERR, name, routine, old_bounds ) 697end if 698if (NEEDS_ALLOC) then 699 allocate( array(b(1,1):b(2,1),b(1,2):b(2,2)), stat=IERR ) 700 call alloc_err( IERR, name, routine, new_bounds ) 701 call alloc_count( size(array), type, name, routine ) 702 array = 0._sp 703end if 704if (NEEDS_COPY) then 705 array(c(1,1):c(2,1),c(1,2):c(2,2)) = & 706 old_array(c(1,1):c(2,1),c(1,2):c(2,2)) 707 call alloc_count( -size(old_array), type, name, routine ) 708 deallocate(old_array,stat=IERR) 709 call alloc_err( IERR, name, routine, old_bounds ) 710end if 711END SUBROUTINE realloc_r2 712! ================================================================== 713SUBROUTINE realloc_r3( array, i1min,i1max, i2min,i2max, i3min,i3max, & 714 name, routine, copy, shrink ) 715implicit none 716character, parameter :: type='R' 717integer, parameter :: rank=3 718real(SP), dimension(:,:,:), pointer :: array, old_array 719integer, intent(in) :: i1min,i1max, i2min,i2max, & 720 i3min,i3max 721character(len=*), optional, intent(in) :: name, routine 722logical, optional, intent(in) :: copy, shrink 723integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 724ASSOCIATED_ARRAY = associated(array) 725if (ASSOCIATED_ARRAY) then 726 old_array => array ! Keep pointer to old array 727 old_bounds(1,:) = lbound(old_array) 728 old_bounds(2,:) = ubound(old_array) 729end if 730new_bounds(1,:) = (/ i1min, i2min, i3min /) 731new_bounds(2,:) = (/ i1max, i2max, i3max /) 732call options( b, c, old_bounds, new_bounds, copy, shrink ) 733if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 734 call alloc_count( -size(old_array), type, name, routine ) 735 deallocate(old_array,stat=IERR) 736 call alloc_err( IERR, name, routine, old_bounds ) 737end if 738if (NEEDS_ALLOC) then 739 allocate( array(b(1,1):b(2,1),b(1,2):b(2,2),b(1,3):b(2,3)),stat=IERR) 740 call alloc_err( IERR, name, routine, new_bounds ) 741 call alloc_count( size(array), type, name, routine ) 742 array = 0._sp 743end if 744if (NEEDS_COPY) then 745 array(c(1,1):c(2,1),c(1,2):c(2,2),c(1,3):c(2,3)) = & 746 old_array(c(1,1):c(2,1),c(1,2):c(2,2),c(1,3):c(2,3)) 747 call alloc_count( -size(old_array), type, name, routine ) 748 deallocate(old_array,stat=IERR) 749 call alloc_err( IERR, name, routine, old_bounds ) 750end if 751END SUBROUTINE realloc_r3 752! ================================================================== 753SUBROUTINE realloc_r4( array, i1min,i1max, i2min,i2max, & 754 i3min,i3max, i4min,i4max, & 755 name, routine, copy, shrink ) 756implicit none 757character, parameter :: type='R' 758integer, parameter :: rank=4 759real(SP), dimension(:,:,:,:), pointer :: array, old_array 760integer, intent(in) :: i1min,i1max, i2min,i2max, & 761 i3min,i3max, i4min,i4max 762character(len=*), optional, intent(in) :: name, routine 763logical, optional, intent(in) :: copy, shrink 764integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 765ASSOCIATED_ARRAY = associated(array) 766if (ASSOCIATED_ARRAY) then 767 old_array => array 768 old_bounds(1,:) = lbound(old_array) 769 old_bounds(2,:) = ubound(old_array) 770end if 771new_bounds(1,:) = (/ i1min, i2min, i3min, i4min /) 772new_bounds(2,:) = (/ i1max, i2max, i3max, i4max /) 773call options( b, c, old_bounds, new_bounds, copy, shrink ) 774if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 775 call alloc_count( -size(old_array), type, name, routine ) 776 deallocate(old_array,stat=IERR) 777 call alloc_err( IERR, name, routine, old_bounds ) 778end if 779if (NEEDS_ALLOC) then 780 allocate( array(b(1,1):b(2,1),b(1,2):b(2,2), & 781 b(1,3):b(2,3),b(1,4):b(2,4)),stat=IERR) 782 call alloc_err( IERR, name, routine, new_bounds ) 783 call alloc_count( size(array), type, name, routine ) 784 array = 0._sp 785end if 786if (NEEDS_COPY) then 787 array(c(1,1):c(2,1),c(1,2):c(2,2),c(1,3):c(2,3),c(1,4):c(2,4))= & 788 old_array(c(1,1):c(2,1),c(1,2):c(2,2),c(1,3):c(2,3),c(1,4):c(2,4)) 789 call alloc_count( -size(old_array), type, name, routine ) 790 deallocate(old_array,stat=IERR) 791 call alloc_err( IERR, name, routine, old_bounds ) 792end if 793END SUBROUTINE realloc_r4 794! ================================================================== 795! Double precision real array reallocs 796! ================================================================== 797SUBROUTINE realloc_d1( array, i1min, i1max, & 798 name, routine, copy, shrink ) 799implicit none 800character, parameter :: type='D' 801integer, parameter :: rank=1 802real(DP), dimension(:), pointer :: array, old_array 803integer, intent(in) :: i1min, i1max 804character(len=*), optional, intent(in) :: name, routine 805logical, optional, intent(in) :: copy, shrink 806integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 807ASSOCIATED_ARRAY = associated(array) 808if (ASSOCIATED_ARRAY) then 809 old_array => array 810 old_bounds(1,:) = lbound(old_array) 811 old_bounds(2,:) = ubound(old_array) 812end if 813new_bounds(1,:) = (/ i1min /) 814new_bounds(2,:) = (/ i1max /) 815call options( b, c, old_bounds, new_bounds, copy, shrink ) 816if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 817 call alloc_count( -size(old_array), type, name, routine ) 818 deallocate(old_array,stat=IERR) 819 call alloc_err( IERR, name, routine, old_bounds ) 820end if 821if (NEEDS_ALLOC) then 822 allocate( array(b(1,1):b(2,1)), stat=IERR ) 823 call alloc_err( IERR, name, routine, new_bounds ) 824 call alloc_count( size(array), type, name, routine ) 825 array = 0._dp 826end if 827if (NEEDS_COPY) then 828 array(c(1,1):c(2,1)) = old_array(c(1,1):c(2,1)) 829 call alloc_count( -size(old_array), type, name, routine ) 830 deallocate(old_array,stat=IERR) 831 call alloc_err( IERR, name, routine, old_bounds ) 832end if 833END SUBROUTINE realloc_d1 834! ================================================================== 835SUBROUTINE realloc_d2( array, i1min,i1max, i2min,i2max, & 836 name, routine, copy, shrink ) 837implicit none 838character, parameter :: type='D' 839integer, parameter :: rank=2 840real(DP), dimension(:,:), pointer :: array, old_array 841integer, intent(in) :: i1min, i1max, i2min, i2max 842character(len=*), optional, intent(in) :: name, routine 843logical, optional, intent(in) :: copy, shrink 844integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 845integer :: i1, i2 846ASSOCIATED_ARRAY = associated(array) 847if (ASSOCIATED_ARRAY) then 848 old_array => array 849 old_bounds(1,:) = lbound(old_array) 850 old_bounds(2,:) = ubound(old_array) 851end if 852new_bounds(1,:) = (/ i1min, i2min /) 853new_bounds(2,:) = (/ i1max, i2max /) 854call options( b, c, old_bounds, new_bounds, copy, shrink ) 855if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 856 call alloc_count( -size(old_array), type, name, routine ) 857 deallocate(old_array,stat=IERR) 858 call alloc_err( IERR, name, routine, old_bounds ) 859end if 860if (NEEDS_ALLOC) then 861 allocate( array(b(1,1):b(2,1),b(1,2):b(2,2)), stat=IERR ) 862 call alloc_err( IERR, name, routine, new_bounds ) 863 call alloc_count( size(array), type, name, routine ) 864 array = 0._dp 865end if 866if (NEEDS_COPY) then 867! array(c(1,1):c(2,1),c(1,2):c(2,2)) = & 868! old_array(c(1,1):c(2,1),c(1,2):c(2,2)) 869 do i2 = c(1,2),c(2,2) 870 do i1 = c(1,1),c(2,1) 871 array(i1,i2) = old_array(i1,i2) 872 end do 873 end do 874 call alloc_count( -size(old_array), type, name, routine ) 875 deallocate(old_array,stat=IERR) 876 call alloc_err( IERR, name, routine, old_bounds ) 877end if 878END SUBROUTINE realloc_d2 879! ================================================================== 880SUBROUTINE realloc_d3( array, i1min,i1max, i2min,i2max, i3min,i3max, & 881 name, routine, copy, shrink ) 882implicit none 883character, parameter :: type='D' 884integer, parameter :: rank=3 885real(DP), dimension(:,:,:), pointer :: array, old_array 886integer, intent(in) :: i1min,i1max, i2min,i2max, & 887 i3min,i3max 888character(len=*), optional, intent(in) :: name, routine 889logical, optional, intent(in) :: copy, shrink 890integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 891integer :: i1, i2, i3 892ASSOCIATED_ARRAY = associated(array) 893if (ASSOCIATED_ARRAY) then 894 old_array => array 895 old_bounds(1,:) = lbound(old_array) 896 old_bounds(2,:) = ubound(old_array) 897end if 898new_bounds(1,:) = (/ i1min, i2min, i3min /) 899new_bounds(2,:) = (/ i1max, i2max, i3max /) 900call options( b, c, old_bounds, new_bounds, copy, shrink ) 901if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 902 call alloc_count( -size(old_array), type, name, routine ) 903 deallocate(old_array,stat=IERR) 904 call alloc_err( IERR, name, routine, old_bounds ) 905end if 906if (NEEDS_ALLOC) then 907 allocate( array(b(1,1):b(2,1),b(1,2):b(2,2),b(1,3):b(2,3)),stat=IERR) 908 call alloc_err( IERR, name, routine, new_bounds ) 909 call alloc_count( size(array), type, name, routine ) 910 array = 0._dp 911end if 912if (NEEDS_COPY) then 913! array(c(1,1):c(2,1),c(1,2):c(2,2),c(1,3):c(2,3)) = & 914! old_array(c(1,1):c(2,1),c(1,2):c(2,2),c(1,3):c(2,3)) 915 do i3 = c(1,3),c(2,3) 916 do i2 = c(1,2),c(2,2) 917 do i1 = c(1,1),c(2,1) 918 array(i1,i2,i3) = old_array(i1,i2,i3) 919 end do 920 end do 921 end do 922 call alloc_count( -size(old_array), type, name, routine ) 923 deallocate(old_array,stat=IERR) 924 call alloc_err( IERR, name, routine, old_bounds ) 925end if 926END SUBROUTINE realloc_d3 927! ================================================================== 928SUBROUTINE realloc_d4( array, i1min,i1max, i2min,i2max, & 929 i3min,i3max, i4min,i4max, & 930 name, routine, copy, shrink ) 931implicit none 932character, parameter :: type='D' 933integer, parameter :: rank=4 934real(DP), dimension(:,:,:,:), pointer :: array, old_array 935integer, intent(in) :: i1min,i1max, i2min,i2max, & 936 i3min,i3max, i4min,i4max 937character(len=*), optional, intent(in) :: name, routine 938logical, optional, intent(in) :: copy, shrink 939integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 940ASSOCIATED_ARRAY = associated(array) 941if (ASSOCIATED_ARRAY) then 942 old_array => array 943 old_bounds(1,:) = lbound(old_array) 944 old_bounds(2,:) = ubound(old_array) 945end if 946new_bounds(1,:) = (/ i1min, i2min, i3min, i4min /) 947new_bounds(2,:) = (/ i1max, i2max, i3max, i4max /) 948call options( b, c, old_bounds, new_bounds, copy, shrink ) 949if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 950 call alloc_count( -size(old_array), type, name, routine ) 951 deallocate(old_array,stat=IERR) 952 call alloc_err( IERR, name, routine, old_bounds ) 953end if 954if (NEEDS_ALLOC) then 955 allocate( array(b(1,1):b(2,1),b(1,2):b(2,2), & 956 b(1,3):b(2,3),b(1,4):b(2,4)),stat=IERR) 957 call alloc_err( IERR, name, routine, new_bounds ) 958 call alloc_count( size(array), type, name, routine ) 959 array = 0._dp 960end if 961if (NEEDS_COPY) then 962 array(c(1,1):c(2,1),c(1,2):c(2,2),c(1,3):c(2,3),c(1,4):c(2,4))= & 963 old_array(c(1,1):c(2,1),c(1,2):c(2,2),c(1,3):c(2,3),c(1,4):c(2,4)) 964 call alloc_count( -size(old_array), type, name, routine ) 965 deallocate(old_array,stat=IERR) 966 call alloc_err( IERR, name, routine, old_bounds ) 967end if 968END SUBROUTINE realloc_d4 969 970! ================================================================== 971! Single precision complex array reallocs 972! ================================================================== 973SUBROUTINE realloc_c1( array, i1min, i1max, & 974 name, routine, copy, shrink ) 975implicit none 976character, parameter :: type='S' 977integer, parameter :: rank=1 978complex(SP), dimension(:), pointer :: array, old_array 979integer, intent(in) :: i1min, i1max 980character(len=*), optional, intent(in) :: name, routine 981logical, optional, intent(in) :: copy, shrink 982integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 983ASSOCIATED_ARRAY = associated(array) 984if (ASSOCIATED_ARRAY) then 985 old_array => array 986 old_bounds(1,:) = lbound(old_array) 987 old_bounds(2,:) = ubound(old_array) 988end if 989new_bounds(1,:) = (/ i1min /) 990new_bounds(2,:) = (/ i1max /) 991call options( b, c, old_bounds, new_bounds, copy, shrink ) 992if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 993 call alloc_count( -2*size(old_array), type, name, routine ) 994 deallocate(old_array,stat=IERR) 995 call alloc_err( IERR, name, routine, old_bounds ) 996end if 997if (NEEDS_ALLOC) then 998 allocate( array(b(1,1):b(2,1)), stat=IERR ) 999 call alloc_err( IERR, name, routine, new_bounds ) 1000 call alloc_count( 2*size(array), type, name, routine ) 1001 array = 0._dp 1002end if 1003if (NEEDS_COPY) then 1004 array(c(1,1):c(2,1)) = old_array(c(1,1):c(2,1)) 1005 call alloc_count( -2*size(old_array), type, name, routine ) 1006 deallocate(old_array,stat=IERR) 1007 call alloc_err( IERR, name, routine, old_bounds ) 1008end if 1009END SUBROUTINE realloc_c1 1010! ================================================================== 1011SUBROUTINE realloc_c2( array, i1min,i1max, i2min,i2max, & 1012 name, routine, copy, shrink ) 1013implicit none 1014character, parameter :: type='S' 1015integer, parameter :: rank=2 1016complex(SP), dimension(:,:), pointer :: array, old_array 1017integer, intent(in) :: i1min, i1max, i2min, i2max 1018character(len=*), optional, intent(in) :: name, routine 1019logical, optional, intent(in) :: copy, shrink 1020integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 1021integer :: i1, i2 1022ASSOCIATED_ARRAY = associated(array) 1023if (ASSOCIATED_ARRAY) then 1024 old_array => array 1025 old_bounds(1,:) = lbound(old_array) 1026 old_bounds(2,:) = ubound(old_array) 1027end if 1028new_bounds(1,:) = (/ i1min, i2min /) 1029new_bounds(2,:) = (/ i1max, i2max /) 1030call options( b, c, old_bounds, new_bounds, copy, shrink ) 1031if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 1032 call alloc_count( -2*size(old_array), type, name, routine ) 1033 deallocate(old_array,stat=IERR) 1034 call alloc_err( IERR, name, routine, old_bounds ) 1035end if 1036if (NEEDS_ALLOC) then 1037 allocate( array(b(1,1):b(2,1),b(1,2):b(2,2)), stat=IERR ) 1038 call alloc_err( IERR, name, routine, new_bounds ) 1039 call alloc_count( 2*size(array), type, name, routine ) 1040 array = 0._dp 1041end if 1042if (NEEDS_COPY) then 1043! array(c(1,1):c(2,1),c(1,2):c(2,2)) = & 1044! old_array(c(1,1):c(2,1),c(1,2):c(2,2)) 1045 do i2 = c(1,2),c(2,2) 1046 do i1 = c(1,1),c(2,1) 1047 array(i1,i2) = old_array(i1,i2) 1048 end do 1049 end do 1050 call alloc_count( -2*size(old_array), type, name, routine ) 1051 deallocate(old_array,stat=IERR) 1052 call alloc_err( IERR, name, routine, old_bounds ) 1053end if 1054END SUBROUTINE realloc_c2 1055 1056! ================================================================== 1057! Double precision complex array reallocs 1058! ================================================================== 1059SUBROUTINE realloc_z1( array, i1min, i1max, & 1060 name, routine, copy, shrink ) 1061implicit none 1062character, parameter :: type='D' 1063integer, parameter :: rank=1 1064complex(DP), dimension(:), pointer :: array, old_array 1065integer, intent(in) :: i1min, i1max 1066character(len=*), optional, intent(in) :: name, routine 1067logical, optional, intent(in) :: copy, shrink 1068integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 1069ASSOCIATED_ARRAY = associated(array) 1070if (ASSOCIATED_ARRAY) then 1071 old_array => array 1072 old_bounds(1,:) = lbound(old_array) 1073 old_bounds(2,:) = ubound(old_array) 1074end if 1075new_bounds(1,:) = (/ i1min /) 1076new_bounds(2,:) = (/ i1max /) 1077call options( b, c, old_bounds, new_bounds, copy, shrink ) 1078if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 1079 call alloc_count( -2*size(old_array), type, name, routine ) 1080 deallocate(old_array,stat=IERR) 1081 call alloc_err( IERR, name, routine, old_bounds ) 1082end if 1083if (NEEDS_ALLOC) then 1084 allocate( array(b(1,1):b(2,1)), stat=IERR ) 1085 call alloc_err( IERR, name, routine, new_bounds ) 1086 call alloc_count( 2*size(array), type, name, routine ) 1087 array = 0._dp 1088end if 1089if (NEEDS_COPY) then 1090 array(c(1,1):c(2,1)) = old_array(c(1,1):c(2,1)) 1091 call alloc_count( -2*size(old_array), type, name, routine ) 1092 deallocate(old_array,stat=IERR) 1093 call alloc_err( IERR, name, routine, old_bounds ) 1094end if 1095END SUBROUTINE realloc_z1 1096! ================================================================== 1097SUBROUTINE realloc_z2( array, i1min,i1max, i2min,i2max, & 1098 name, routine, copy, shrink ) 1099implicit none 1100character, parameter :: type='D' 1101integer, parameter :: rank=2 1102complex(DP), dimension(:,:), pointer :: array, old_array 1103integer, intent(in) :: i1min, i1max, i2min, i2max 1104character(len=*), optional, intent(in) :: name, routine 1105logical, optional, intent(in) :: copy, shrink 1106integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 1107integer :: i1, i2 1108ASSOCIATED_ARRAY = associated(array) 1109if (ASSOCIATED_ARRAY) then 1110 old_array => array 1111 old_bounds(1,:) = lbound(old_array) 1112 old_bounds(2,:) = ubound(old_array) 1113end if 1114new_bounds(1,:) = (/ i1min, i2min /) 1115new_bounds(2,:) = (/ i1max, i2max /) 1116call options( b, c, old_bounds, new_bounds, copy, shrink ) 1117if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 1118 call alloc_count( -2*size(old_array), type, name, routine ) 1119 deallocate(old_array,stat=IERR) 1120 call alloc_err( IERR, name, routine, old_bounds ) 1121end if 1122if (NEEDS_ALLOC) then 1123 allocate( array(b(1,1):b(2,1),b(1,2):b(2,2)), stat=IERR ) 1124 call alloc_err( IERR, name, routine, new_bounds ) 1125 call alloc_count( 2*size(array), type, name, routine ) 1126 array = 0._dp 1127end if 1128if (NEEDS_COPY) then 1129! array(c(1,1):c(2,1),c(1,2):c(2,2)) = & 1130! old_array(c(1,1):c(2,1),c(1,2):c(2,2)) 1131 do i2 = c(1,2),c(2,2) 1132 do i1 = c(1,1),c(2,1) 1133 array(i1,i2) = old_array(i1,i2) 1134 end do 1135 end do 1136 call alloc_count( -2*size(old_array), type, name, routine ) 1137 deallocate(old_array,stat=IERR) 1138 call alloc_err( IERR, name, routine, old_bounds ) 1139end if 1140END SUBROUTINE realloc_z2 1141SUBROUTINE realloc_z3( array, i1min,i1max, i2min,i2max, i3min,i3max, & 1142 name, routine, copy, shrink ) 1143implicit none 1144character, parameter :: type='D' 1145integer, parameter :: rank=3 1146complex(DP), dimension(:,:,:), pointer :: array, old_array 1147integer, intent(in) :: i1min, i1max, i2min, i2max 1148integer, intent(in) :: i3min, i3max 1149character(len=*), optional, intent(in) :: name, routine 1150logical, optional, intent(in) :: copy, shrink 1151integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 1152integer :: i1, i2, i3 1153ASSOCIATED_ARRAY = associated(array) 1154if (ASSOCIATED_ARRAY) then 1155 old_array => array 1156 old_bounds(1,:) = lbound(old_array) 1157 old_bounds(2,:) = ubound(old_array) 1158end if 1159new_bounds(1,:) = (/ i1min, i2min, i3min /) 1160new_bounds(2,:) = (/ i1max, i2max, i3max /) 1161call options( b, c, old_bounds, new_bounds, copy, shrink ) 1162if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 1163 call alloc_count( -2*size(old_array), type, name, routine ) 1164 deallocate(old_array,stat=IERR) 1165 call alloc_err( IERR, name, routine, old_bounds ) 1166end if 1167if (NEEDS_ALLOC) then 1168 allocate( array(b(1,1):b(2,1),b(1,2):b(2,2),b(1,3):b(2,3)), stat=IERR ) 1169 call alloc_err( IERR, name, routine, new_bounds ) 1170 call alloc_count( 2*size(array), type, name, routine ) 1171 array = 0._dp 1172end if 1173if (NEEDS_COPY) then 1174! array(c(1,1):c(2,1),c(1,2):c(2,2)) = & 1175! old_array(c(1,1):c(2,1),c(1,2):c(2,2)) 1176 do i3 = c(1,3),c(2,3) 1177 do i2 = c(1,2),c(2,2) 1178 do i1 = c(1,1),c(2,1) 1179 array(i1,i2,i3) = old_array(i1,i2,i3) 1180 end do 1181 end do 1182 end do 1183 call alloc_count( -2*size(old_array), type, name, routine ) 1184 deallocate(old_array,stat=IERR) 1185 call alloc_err( IERR, name, routine, old_bounds ) 1186end if 1187END SUBROUTINE realloc_z3 1188SUBROUTINE realloc_z4( array, i1min,i1max, i2min,i2max, & 1189 i3min,i3max, i4min,i4max, & 1190 name, routine, copy, shrink ) 1191implicit none 1192character, parameter :: type='D' 1193integer, parameter :: rank=4 1194complex(DP), dimension(:,:,:,:), pointer :: array, old_array 1195integer, intent(in) :: i1min, i1max, i2min, i2max, & 1196 i3min, i3max, i4min, i4max 1197character(len=*), optional, intent(in) :: name, routine 1198logical, optional, intent(in) :: copy, shrink 1199integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 1200integer :: i1, i2, i3, i4 1201ASSOCIATED_ARRAY = associated(array) 1202if (ASSOCIATED_ARRAY) then 1203 old_array => array 1204 old_bounds(1,:) = lbound(old_array) 1205 old_bounds(2,:) = ubound(old_array) 1206end if 1207new_bounds(1,:) = (/ i1min, i2min, i3min, i4min /) 1208new_bounds(2,:) = (/ i1max, i2max, i3max, i4max /) 1209call options( b, c, old_bounds, new_bounds, copy, shrink ) 1210if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 1211 call alloc_count( -2*size(old_array), type, name, routine ) 1212 deallocate(old_array,stat=IERR) 1213 call alloc_err( IERR, name, routine, old_bounds ) 1214end if 1215if (NEEDS_ALLOC) then 1216 allocate( array(b(1,1):b(2,1),b(1,2):b(2,2),b(1,3):b(2,3),b(1,4):b(2,4)), & 1217 stat=IERR ) 1218 call alloc_err( IERR, name, routine, new_bounds ) 1219 call alloc_count( 2*size(array), type, name, routine ) 1220 array = 0._dp 1221end if 1222if (NEEDS_COPY) then 1223! array(c(1,1):c(2,1),c(1,2):c(2,2),c(1,3):c(2,3),c(1,4):c(2,4)) = & 1224! old_array(c(1,1):c(2,1),c(1,2):c(2,2),c(1,3):c(2,3),c(1,4):c(2,4)) 1225 do i4 = c(1,4),c(2,4) 1226 do i3 = c(1,3),c(2,3) 1227 do i2 = c(1,2),c(2,2) 1228 do i1 = c(1,1),c(2,1) 1229 array(i1,i2,i3,i4) = old_array(i1,i2,i3,i4) 1230 end do 1231 end do 1232 end do 1233 end do 1234 call alloc_count( -2*size(old_array), type, name, routine ) 1235 deallocate(old_array,stat=IERR) 1236 call alloc_err( IERR, name, routine, old_bounds ) 1237end if 1238END SUBROUTINE realloc_z4 1239! ================================================================== 1240! Logical array reallocs 1241! ================================================================== 1242SUBROUTINE realloc_l1( array, i1min,i1max, & 1243 name, routine, copy, shrink ) 1244implicit none 1245character, parameter :: type='L' 1246integer, parameter :: rank=1 1247logical, dimension(:), pointer :: array, old_array 1248integer, intent(in) :: i1min,i1max 1249character(len=*), optional, intent(in) :: name, routine 1250logical, optional, intent(in) :: copy, shrink 1251integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 1252ASSOCIATED_ARRAY = associated(array) 1253if (ASSOCIATED_ARRAY) then 1254 old_array => array 1255 old_bounds(1,:) = lbound(old_array) 1256 old_bounds(2,:) = ubound(old_array) 1257end if 1258new_bounds(1,:) = (/ i1min /) 1259new_bounds(2,:) = (/ i1max /) 1260call options( b, c, old_bounds, new_bounds, copy, shrink ) 1261if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 1262 call alloc_count( -size(old_array), type, name, routine ) 1263 deallocate(old_array,stat=IERR) 1264 call alloc_err( IERR, name, routine, old_bounds ) 1265end if 1266if (NEEDS_ALLOC) then 1267 allocate( array(b(1,1):b(2,1)), stat=IERR ) 1268 call alloc_err( IERR, name, routine, new_bounds ) 1269 call alloc_count( size(array), type, name, routine ) 1270 array = .false. 1271end if 1272if (NEEDS_COPY) then 1273 array(c(1,1):c(2,1)) = old_array(c(1,1):c(2,1)) 1274 call alloc_count( -size(old_array), type, name, routine ) 1275 deallocate(old_array,stat=IERR) 1276 call alloc_err( IERR, name, routine, old_bounds ) 1277end if 1278END SUBROUTINE realloc_l1 1279! ================================================================== 1280SUBROUTINE realloc_l2( array, i1min,i1max, i2min,i2max, & 1281 name, routine, copy, shrink ) 1282implicit none 1283character, parameter :: type='L' 1284integer, parameter :: rank=2 1285logical, dimension(:,:), pointer :: array, old_array 1286integer, intent(in) :: i1min,i1max, i2min,i2max 1287character(len=*), optional, intent(in) :: name, routine 1288logical, optional, intent(in) :: copy, shrink 1289integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 1290ASSOCIATED_ARRAY = associated(array) 1291if (ASSOCIATED_ARRAY) then 1292 old_array => array 1293 old_bounds(1,:) = lbound(old_array) 1294 old_bounds(2,:) = ubound(old_array) 1295end if 1296new_bounds(1,:) = (/ i1min, i2min /) 1297new_bounds(2,:) = (/ i1max, i2max /) 1298call options( b, c, old_bounds, new_bounds, copy, shrink ) 1299if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 1300 call alloc_count( -size(old_array), type, name, routine ) 1301 deallocate(old_array,stat=IERR) 1302 call alloc_err( IERR, name, routine, old_bounds ) 1303end if 1304if (NEEDS_ALLOC) then 1305 allocate( array(b(1,1):b(2,1),b(1,2):b(2,2)), stat=IERR ) 1306 call alloc_err( IERR, name, routine, new_bounds ) 1307 call alloc_count( size(array), type, name, routine ) 1308 array = .false. 1309end if 1310if (NEEDS_COPY) then 1311 array(c(1,1):c(2,1),c(1,2):c(2,2)) = & 1312 old_array(c(1,1):c(2,1),c(1,2):c(2,2)) 1313 call alloc_count( -size(old_array), type, name, routine ) 1314 deallocate(old_array,stat=IERR) 1315 call alloc_err( IERR, name, routine, old_bounds ) 1316end if 1317END SUBROUTINE realloc_l2 1318! ================================================================== 1319SUBROUTINE realloc_l3( array, i1min,i1max, i2min,i2max, i3min,i3max, & 1320 name, routine, copy, shrink ) 1321implicit none 1322character, parameter :: type='L' 1323integer, parameter :: rank=3 1324logical, dimension(:,:,:), pointer :: array, old_array 1325integer, intent(in) :: i1min,i1max, i2min,i2max, & 1326 i3min,i3max 1327character(len=*), optional, intent(in) :: name, routine 1328logical, optional, intent(in) :: copy, shrink 1329integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 1330ASSOCIATED_ARRAY = associated(array) 1331if (ASSOCIATED_ARRAY) then 1332 old_array => array 1333 old_bounds(1,:) = lbound(old_array) 1334 old_bounds(2,:) = ubound(old_array) 1335end if 1336new_bounds(1,:) = (/ i1min, i2min, i3min /) 1337new_bounds(2,:) = (/ i1max, i2max, i3max /) 1338call options( b, c, old_bounds, new_bounds, copy, shrink ) 1339if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 1340 call alloc_count( -size(old_array), type, name, routine ) 1341 deallocate(old_array,stat=IERR) 1342 call alloc_err( IERR, name, routine, old_bounds ) 1343end if 1344if (NEEDS_ALLOC) then 1345 allocate( array(b(1,1):b(2,1),b(1,2):b(2,2),b(1,3):b(2,3)),stat=IERR) 1346 call alloc_err( IERR, name, routine, new_bounds ) 1347 call alloc_count( size(array), type, name, routine ) 1348 array = .false. 1349end if 1350if (NEEDS_COPY) then 1351 array(c(1,1):c(2,1),c(1,2):c(2,2),c(1,3):c(2,3)) = & 1352 old_array(c(1,1):c(2,1),c(1,2):c(2,2),c(1,3):c(2,3)) 1353 call alloc_count( -size(old_array), type, name, routine ) 1354 deallocate(old_array,stat=IERR) 1355 call alloc_err( IERR, name, routine, old_bounds ) 1356end if 1357END SUBROUTINE realloc_l3 1358! ================================================================== 1359! Realloc routines with assumed lower bound = 1 1360!AG: Extremely dangerous -- do not use. 1361! ================================================================== 1362SUBROUTINE realloc_i1s( array, i1max, & 1363 name, routine, copy, shrink ) 1364! Arguments 1365implicit none 1366integer, dimension(:), pointer :: array 1367integer, intent(in) :: i1max 1368character(len=*), optional, intent(in) :: name 1369character(len=*), optional, intent(in) :: routine 1370logical, optional, intent(in) :: copy 1371logical, optional, intent(in) :: shrink 1372 1373call realloc_i1( array, DEFAULT%imin, i1max, & 1374 name, routine, copy, shrink ) 1375 1376END SUBROUTINE realloc_i1s 1377! ================================================================== 1378SUBROUTINE realloc_i2s( array, i1max, i2max, & 1379 name, routine, copy, shrink ) 1380implicit none 1381integer, dimension(:,:), pointer :: array 1382integer, intent(in) :: i1max, i2max 1383character(len=*), optional, intent(in) :: name, routine 1384logical, optional, intent(in) :: copy, shrink 1385call realloc_i2( array, DEFAULT%imin, i1max, DEFAULT%imin, i2max, & 1386 name, routine, copy, shrink ) 1387END SUBROUTINE realloc_i2s 1388! ================================================================== 1389SUBROUTINE realloc_i3s( array, i1max, i2max, i3max, & 1390 name, routine, copy, shrink ) 1391implicit none 1392integer, dimension(:,:,:), pointer :: array 1393integer, intent(in) :: i1max, i2max, i3max 1394character(len=*), optional, intent(in) :: name, routine 1395logical, optional, intent(in) :: copy, shrink 1396call realloc_i3( array, DEFAULT%imin, i1max, DEFAULT%imin, i2max, & 1397 DEFAULT%imin, i3max, & 1398 name, routine, copy, shrink ) 1399END SUBROUTINE realloc_i3s 1400! ================================================================== 1401SUBROUTINE realloc_r1s( array, i1max, & 1402 name, routine, copy, shrink ) 1403implicit none 1404real(SP), dimension(:), pointer :: array 1405integer, intent(in) :: i1max 1406character(len=*), optional, intent(in) :: name, routine 1407logical, optional, intent(in) :: copy, shrink 1408call realloc_r1( array, DEFAULT%imin, i1max, & 1409 name, routine, copy, shrink ) 1410END SUBROUTINE realloc_r1s 1411! ================================================================== 1412SUBROUTINE realloc_r2s( array, i1max, i2max, & 1413 name, routine, copy, shrink ) 1414implicit none 1415real(SP), dimension(:,:), pointer :: array 1416integer, intent(in) :: i1max, i2max 1417character(len=*), optional, intent(in) :: name, routine 1418logical, optional, intent(in) :: copy, shrink 1419call realloc_r2( array, DEFAULT%imin, i1max, DEFAULT%imin, i2max, & 1420 name, routine, copy, shrink ) 1421END SUBROUTINE realloc_r2s 1422! ================================================================== 1423SUBROUTINE realloc_r3s( array, i1max, i2max, i3max, & 1424 name, routine, copy, shrink ) 1425implicit none 1426real(SP), dimension(:,:,:), pointer :: array 1427integer, intent(in) :: i1max, i2max, i3max 1428character(len=*), optional, intent(in) :: name, routine 1429logical, optional, intent(in) :: copy, shrink 1430call realloc_r3( array, DEFAULT%imin, i1max, DEFAULT%imin, i2max, & 1431 DEFAULT%imin, i3max, & 1432 name, routine, copy, shrink ) 1433END SUBROUTINE realloc_r3s 1434! ================================================================== 1435SUBROUTINE realloc_r4s( array, i1max, i2max, i3max, i4max, & 1436 name, routine, copy, shrink ) 1437implicit none 1438real(SP), dimension(:,:,:,:), pointer :: array 1439integer, intent(in) :: i1max, i2max, i3max, i4max 1440character(len=*), optional, intent(in) :: name, routine 1441logical, optional, intent(in) :: copy, shrink 1442call realloc_r4( array, DEFAULT%imin, i1max, DEFAULT%imin, i2max, & 1443 DEFAULT%imin, i3max, DEFAULT%imin, i4max, & 1444 name, routine, copy, shrink ) 1445END SUBROUTINE realloc_r4s 1446! ================================================================== 1447SUBROUTINE realloc_d1s( array, i1max, & 1448 name, routine, copy, shrink ) 1449implicit none 1450real(DP), dimension(:), pointer :: array 1451integer, intent(in) :: i1max 1452character(len=*), optional, intent(in) :: name, routine 1453logical, optional, intent(in) :: copy, shrink 1454call realloc_d1( array, DEFAULT%imin, i1max, & 1455 name, routine, copy, shrink ) 1456END SUBROUTINE realloc_d1s 1457! ================================================================== 1458SUBROUTINE realloc_d2s( array, i1max, i2max, & 1459 name, routine, copy, shrink ) 1460implicit none 1461real(DP), dimension(:,:), pointer :: array 1462integer, intent(in) :: i1max, i2max 1463character(len=*), optional, intent(in) :: name, routine 1464logical, optional, intent(in) :: copy, shrink 1465call realloc_d2( array, DEFAULT%imin, i1max, DEFAULT%imin, i2max, & 1466 name, routine, copy, shrink ) 1467END SUBROUTINE realloc_d2s 1468! ================================================================== 1469SUBROUTINE realloc_d3s( array, i1max, i2max, i3max, & 1470 name, routine, copy, shrink ) 1471implicit none 1472real(DP), dimension(:,:,:), pointer :: array 1473integer, intent(in) :: i1max, i2max, i3max 1474character(len=*), optional, intent(in) :: name, routine 1475logical, optional, intent(in) :: copy, shrink 1476call realloc_d3( array, DEFAULT%imin, i1max, DEFAULT%imin, i2max, & 1477 DEFAULT%imin, i3max, & 1478 name, routine, copy, shrink ) 1479END SUBROUTINE realloc_d3s 1480! ================================================================== 1481SUBROUTINE realloc_d4s( array, i1max, i2max, i3max, i4max, & 1482 name, routine, copy, shrink ) 1483implicit none 1484real(DP), dimension(:,:,:,:), pointer :: array 1485integer, intent(in) :: i1max, i2max, i3max, i4max 1486character(len=*), optional, intent(in) :: name, routine 1487logical, optional, intent(in) :: copy, shrink 1488call realloc_d4( array, DEFAULT%imin, i1max, DEFAULT%imin, i2max, & 1489 DEFAULT%imin, i3max, DEFAULT%imin, i4max, & 1490 name, routine, copy, shrink ) 1491END SUBROUTINE realloc_d4s 1492! ================================================================== 1493SUBROUTINE realloc_l1s( array, i1max, & 1494 name, routine, copy, shrink ) 1495implicit none 1496logical, dimension(:), pointer :: array 1497integer, intent(in) :: i1max 1498character(len=*), optional, intent(in) :: name 1499character(len=*), optional, intent(in) :: routine 1500logical, optional, intent(in) :: copy 1501logical, optional, intent(in) :: shrink 1502call realloc_l1( array, DEFAULT%imin, i1max, & 1503 name, routine, copy, shrink ) 1504END SUBROUTINE realloc_l1s 1505! ================================================================== 1506SUBROUTINE realloc_l2s( array, i1max, i2max, & 1507 name, routine, copy, shrink ) 1508implicit none 1509logical, dimension(:,:), pointer :: array 1510integer, intent(in) :: i1max, i2max 1511character(len=*), optional, intent(in) :: name 1512character(len=*), optional, intent(in) :: routine 1513logical, optional, intent(in) :: copy 1514logical, optional, intent(in) :: shrink 1515call realloc_l2( array, DEFAULT%imin, i1max, DEFAULT%imin, i2max, & 1516 name, routine, copy, shrink ) 1517END SUBROUTINE realloc_l2s 1518! ================================================================== 1519SUBROUTINE realloc_l3s( array, i1max, i2max, i3max, & 1520 name, routine, copy, shrink ) 1521implicit none 1522logical, dimension(:,:,:), pointer :: array 1523integer, intent(in) :: i1max, i2max, i3max 1524character(len=*), optional, intent(in) :: name 1525character(len=*), optional, intent(in) :: routine 1526logical, optional, intent(in) :: copy 1527logical, optional, intent(in) :: shrink 1528call realloc_l3( array, DEFAULT%imin, i1max, DEFAULT%imin, i2max, & 1529 DEFAULT%imin, i3max, name, routine, copy, shrink ) 1530END SUBROUTINE realloc_l3s 1531! ================================================================== 1532! Character vector realloc 1533! ================================================================== 1534SUBROUTINE realloc_s1( array, i1min, i1max, & 1535 name, routine, copy, shrink ) 1536! Arguments 1537implicit none 1538character(len=*), dimension(:), pointer :: array 1539integer, intent(in) :: i1min 1540integer, intent(in) :: i1max 1541character(len=*), optional, intent(in) :: name 1542character(len=*), optional, intent(in) :: routine 1543logical, optional, intent(in) :: copy 1544logical, optional, intent(in) :: shrink 1545 1546! Internal variables and arrays 1547character, parameter :: type='S' 1548integer, parameter :: rank=1 1549character(len=len(array)), dimension(:), pointer :: old_array 1550integer, dimension(2,rank) :: b, c, new_bounds, old_bounds 1551 1552! Get old array bounds 1553ASSOCIATED_ARRAY = associated(array) 1554if (ASSOCIATED_ARRAY) then 1555 old_array => array ! Keep pointer to old array 1556 old_bounds(1,:) = lbound(old_array) 1557 old_bounds(2,:) = ubound(old_array) 1558end if 1559 1560! Copy new requested array bounds 1561new_bounds(1,:) = (/ i1min /) 1562new_bounds(2,:) = (/ i1max /) 1563 1564! Find if it is a new allocation or a true reallocation, 1565! and if the contents need to be copied (saved) 1566! Argument b returns common bounds 1567! Options routine also reads common variable ASSOCIATED_ARRAY, 1568! and it sets NEEDS_ALLOC, NEEDS_DEALLOC, and NEEDS_COPY 1569call options( b, c, old_bounds, new_bounds, copy, shrink ) 1570 1571! Deallocate old space 1572if (NEEDS_DEALLOC .and. .not.NEEDS_COPY) then 1573 call alloc_count( -size(old_array)*len(old_array), type, name, routine ) 1574 deallocate(old_array,stat=IERR) 1575 call alloc_err( IERR, name, routine, old_bounds ) 1576end if 1577 1578! Allocate new space 1579if (NEEDS_ALLOC) then 1580 allocate( array(b(1,1):b(2,1)), stat=IERR ) 1581 call alloc_err( IERR, name, routine, new_bounds ) 1582 call alloc_count( size(array)*len(array), type, name, routine ) 1583 array = '' 1584end if 1585 1586! Copy contents and deallocate old space 1587if (NEEDS_COPY) then 1588 array(c(1,1):c(2,1)) = old_array(c(1,1):c(2,1)) 1589 call alloc_count( -size(old_array)*len(old_array), type, name, routine ) 1590 deallocate(old_array,stat=IERR) 1591 call alloc_err( IERR, name, routine, old_bounds ) 1592end if 1593 1594END SUBROUTINE realloc_s1 1595! ================================================================== 1596! Dealloc routines 1597! ================================================================== 1598SUBROUTINE dealloc_i1( array, name, routine ) 1599 1600! Arguments 1601implicit none 1602integer, dimension(:), pointer :: array 1603character(len=*), optional, intent(in) :: name 1604character(len=*), optional, intent(in) :: routine 1605 1606if (associated(array)) then 1607 call alloc_count( -size(array), 'I', name, routine ) 1608 deallocate(array,stat=IERR) 1609 call alloc_err( IERR, name, routine ) 1610end if 1611 1612END SUBROUTINE dealloc_i1 1613 1614! ================================================================== 1615SUBROUTINE dealloc_i2( array, name, routine ) 1616implicit none 1617integer, dimension(:,:), pointer :: array 1618character(len=*), optional, intent(in) :: name, routine 1619if (associated(array)) then 1620 call alloc_count( -size(array), 'I', name, routine ) 1621 deallocate(array,stat=IERR) 1622 call alloc_err( IERR, name, routine ) 1623 1624end if 1625END SUBROUTINE dealloc_i2 1626! ================================================================== 1627SUBROUTINE dealloc_i3( array, name, routine ) 1628implicit none 1629integer, dimension(:,:,:), pointer :: array 1630character(len=*), optional, intent(in) :: name, routine 1631if (associated(array)) then 1632 call alloc_count( -size(array), 'I', name, routine ) 1633 deallocate(array,stat=IERR) 1634 call alloc_err( IERR, name, routine ) 1635 1636end if 1637END SUBROUTINE dealloc_i3 1638! ================================================================== 1639SUBROUTINE dealloc_E1( array, name, routine ) 1640 1641! Arguments 1642implicit none 1643integer(i8b), dimension(:), pointer :: array 1644character(len=*), optional, intent(in) :: name 1645character(len=*), optional, intent(in) :: routine 1646 1647if (associated(array)) then 1648 call alloc_count( -size(array), 'I', name, routine ) 1649 deallocate(array,stat=IERR) 1650 call alloc_err( IERR, name, routine ) 1651 1652end if 1653 1654END SUBROUTINE dealloc_E1 1655 1656! ================================================================== 1657SUBROUTINE dealloc_r1( array, name, routine ) 1658implicit none 1659real(SP), dimension(:), pointer :: array 1660character(len=*), optional, intent(in) :: name, routine 1661if (associated(array)) then 1662 call alloc_count( -size(array), 'R', name, routine ) 1663 deallocate(array,stat=IERR) 1664 call alloc_err( IERR, name, routine ) 1665end if 1666END SUBROUTINE dealloc_r1 1667! ================================================================== 1668SUBROUTINE dealloc_r2( array, name, routine ) 1669implicit none 1670real(SP), dimension(:,:), pointer :: array 1671character(len=*), optional, intent(in) :: name, routine 1672if (associated(array)) then 1673 call alloc_count( -size(array), 'R', name, routine ) 1674 deallocate(array,stat=IERR) 1675 call alloc_err( IERR, name, routine ) 1676end if 1677END SUBROUTINE dealloc_r2 1678! ================================================================== 1679SUBROUTINE dealloc_r3( array, name, routine ) 1680implicit none 1681real(SP), dimension(:,:,:), pointer :: array 1682character(len=*), optional, intent(in) :: name, routine 1683if (associated(array)) then 1684 call alloc_count( -size(array), 'R', name, routine ) 1685 deallocate(array,stat=IERR) 1686 call alloc_err( IERR, name, routine ) 1687end if 1688END SUBROUTINE dealloc_r3 1689! ================================================================== 1690SUBROUTINE dealloc_r4( array, name, routine ) 1691implicit none 1692real(SP), dimension(:,:,:,:), pointer :: array 1693character(len=*), optional, intent(in) :: name, routine 1694if (associated(array)) then 1695 call alloc_count( -size(array), 'R', name, routine ) 1696 deallocate(array,stat=IERR) 1697 call alloc_err( IERR, name, routine ) 1698end if 1699END SUBROUTINE dealloc_r4 1700! ================================================================== 1701SUBROUTINE dealloc_d1( array, name, routine ) 1702implicit none 1703real(DP), dimension(:), pointer :: array 1704character(len=*), optional, intent(in) :: name, routine 1705if (associated(array)) then 1706 call alloc_count( -size(array), 'D', name, routine ) 1707 deallocate(array,stat=IERR) 1708 call alloc_err( IERR, name, routine ) 1709end if 1710END SUBROUTINE dealloc_d1 1711! ================================================================== 1712SUBROUTINE dealloc_d2( array, name, routine ) 1713implicit none 1714real(DP), dimension(:,:), pointer :: array 1715character(len=*), optional, intent(in) :: name, routine 1716if (associated(array)) then 1717 call alloc_count( -size(array), 'D', name, routine ) 1718 deallocate(array,stat=IERR) 1719 call alloc_err( IERR, name, routine ) 1720end if 1721END SUBROUTINE dealloc_d2 1722! ================================================================== 1723SUBROUTINE dealloc_d3( array, name, routine ) 1724implicit none 1725real(DP), dimension(:,:,:), pointer :: array 1726character(len=*), optional, intent(in) :: name, routine 1727if (associated(array)) then 1728 call alloc_count( -size(array), 'D', name, routine ) 1729 deallocate(array,stat=IERR) 1730 call alloc_err( IERR, name, routine ) 1731end if 1732END SUBROUTINE dealloc_d3 1733! ================================================================== 1734SUBROUTINE dealloc_d4( array, name, routine ) 1735implicit none 1736real(DP), dimension(:,:,:,:), pointer :: array 1737character(len=*), optional, intent(in) :: name, routine 1738if (associated(array)) then 1739 call alloc_count( -size(array), 'D', name, routine ) 1740 deallocate(array,stat=IERR) 1741 call alloc_err( IERR, name, routine ) 1742end if 1743END SUBROUTINE dealloc_d4 1744! ================================================================== 1745! COMPLEX versions 1746! 1747SUBROUTINE dealloc_c1( array, name, routine ) 1748implicit none 1749complex(SP), dimension(:), pointer :: array 1750character(len=*), optional, intent(in) :: name, routine 1751if (associated(array)) then 1752 call alloc_count( -2*size(array), 'S', name, routine ) 1753 deallocate(array,stat=IERR) 1754 call alloc_err( IERR, name, routine ) 1755end if 1756END SUBROUTINE dealloc_c1 1757! ================================================================== 1758SUBROUTINE dealloc_c2( array, name, routine ) 1759implicit none 1760complex(SP), dimension(:,:), pointer :: array 1761character(len=*), optional, intent(in) :: name, routine 1762if (associated(array)) then 1763 call alloc_count( -2*size(array), 'S', name, routine ) 1764 deallocate(array,stat=IERR) 1765 call alloc_err( IERR, name, routine ) 1766end if 1767END SUBROUTINE dealloc_c2 1768! ================================================================== 1769SUBROUTINE dealloc_z1( array, name, routine ) 1770implicit none 1771complex(DP), dimension(:), pointer :: array 1772character(len=*), optional, intent(in) :: name, routine 1773if (associated(array)) then 1774 call alloc_count( -2*size(array), 'D', name, routine ) 1775 deallocate(array,stat=IERR) 1776 call alloc_err( IERR, name, routine ) 1777end if 1778END SUBROUTINE dealloc_z1 1779! ================================================================== 1780SUBROUTINE dealloc_z2( array, name, routine ) 1781implicit none 1782complex(DP), dimension(:,:), pointer :: array 1783character(len=*), optional, intent(in) :: name, routine 1784if (associated(array)) then 1785 call alloc_count( -2*size(array), 'D', name, routine ) 1786 deallocate(array,stat=IERR) 1787 call alloc_err( IERR, name, routine ) 1788end if 1789END SUBROUTINE dealloc_z2 1790! ================================================================== 1791SUBROUTINE dealloc_z3( array, name, routine ) 1792implicit none 1793complex(DP), dimension(:,:,:), pointer :: array 1794character(len=*), optional, intent(in) :: name, routine 1795if (associated(array)) then 1796 call alloc_count( -2*size(array), 'D', name, routine ) 1797 deallocate(array,stat=IERR) 1798 call alloc_err( IERR, name, routine ) 1799end if 1800END SUBROUTINE dealloc_z3 1801! ================================================================== 1802SUBROUTINE dealloc_z4( array, name, routine ) 1803implicit none 1804complex(DP), dimension(:,:,:,:), pointer :: array 1805character(len=*), optional, intent(in) :: name, routine 1806if (associated(array)) then 1807 call alloc_count( -2*size(array), 'D', name, routine ) 1808 deallocate(array,stat=IERR) 1809 call alloc_err( IERR, name, routine ) 1810end if 1811END SUBROUTINE dealloc_z4 1812! ================================================================== 1813SUBROUTINE dealloc_l1( array, name, routine ) 1814implicit none 1815logical, dimension(:), pointer :: array 1816character(len=*), optional, intent(in) :: name, routine 1817if (associated(array)) then 1818 call alloc_count( -size(array), 'L', name, routine ) 1819 deallocate(array,stat=IERR) 1820 call alloc_err( IERR, name, routine ) 1821end if 1822END SUBROUTINE dealloc_l1 1823! ================================================================== 1824SUBROUTINE dealloc_l2( array, name, routine ) 1825implicit none 1826logical, dimension(:,:), pointer :: array 1827character(len=*), optional, intent(in) :: name, routine 1828if (associated(array)) then 1829 call alloc_count( -size(array), 'L', name, routine ) 1830 deallocate(array,stat=IERR) 1831 call alloc_err( IERR, name, routine ) 1832end if 1833END SUBROUTINE dealloc_l2 1834! ================================================================== 1835SUBROUTINE dealloc_l3( array, name, routine ) 1836implicit none 1837logical, dimension(:,:,:), pointer :: array 1838character(len=*), optional, intent(in) :: name, routine 1839if (associated(array)) then 1840 call alloc_count( -size(array), 'L', name, routine ) 1841 deallocate(array,stat=IERR) 1842 call alloc_err( IERR, name, routine ) 1843end if 1844END SUBROUTINE dealloc_l3 1845! ================================================================== 1846SUBROUTINE dealloc_s1( array, name, routine ) 1847implicit none 1848character(len=*), dimension(:), pointer :: array 1849character(len=*), optional, intent(in) :: name, routine 1850if (associated(array)) then 1851 call alloc_count( -size(array)*len(array), 'S', name, routine ) 1852 deallocate(array,stat=IERR) 1853 call alloc_err( IERR, name, routine ) 1854end if 1855END SUBROUTINE dealloc_s1 1856 1857! ================================================================== 1858! Internal subroutines 1859! ================================================================== 1860 1861SUBROUTINE options( final_bounds, common_bounds, & 1862 old_bounds, new_bounds, copy, shrink ) 1863! Arguments 1864integer, dimension(:,:), intent(out) :: final_bounds 1865integer, dimension(:,:), intent(out) :: common_bounds 1866integer, dimension(:,:), intent(in) :: old_bounds 1867integer, dimension(:,:), intent(in) :: new_bounds 1868logical, optional, intent(in) :: copy 1869logical, optional, intent(in) :: shrink 1870 1871! Internal variables and arrays 1872logical want_shrink 1873 1874 1875!! AG***** 1876! It might be worthwhile to check whether the user 1877! atttemps to use bounds which do not make sense, 1878! such as zero, or with upper<lower... 1879!!*** 1880 1881! Find if it is a new allocation or a true reallocation, 1882! and if the contents need to be copied (saved) 1883if (ASSOCIATED_ARRAY) then 1884 1885 ! Check if array bounds have changed 1886 if ( all(new_bounds==old_bounds) ) then 1887 ! Old and new arrays are equal. Nothing needs to be done 1888 NEEDS_ALLOC = .false. 1889 NEEDS_DEALLOC = .false. 1890 NEEDS_COPY = .false. 1891 else 1892 1893 ! Want to shrink? 1894 if (present(shrink)) then 1895 want_shrink = shrink 1896 else 1897 want_shrink = DEFAULT%shrink 1898 end if 1899 1900 if (.not. want_shrink & 1901 .and. all(new_bounds(1,:)>=old_bounds(1,:)) & 1902 .and. all(new_bounds(2,:)<=old_bounds(2,:)) ) then 1903 ! Old array is already fine. Nothing needs to be done 1904 NEEDS_ALLOC = .false. 1905 NEEDS_DEALLOC = .false. 1906 NEEDS_COPY = .false. 1907 else 1908 ! Old array needs to be substituted by a new array 1909 NEEDS_ALLOC = .true. 1910 NEEDS_DEALLOC = .true. 1911 if (present(copy)) then 1912 NEEDS_COPY = copy 1913 else 1914 NEEDS_COPY = DEFAULT%copy 1915 end if 1916 1917 ! Ensure that bounds shrink only if desired 1918 if (want_shrink) then 1919 final_bounds(1,:) = new_bounds(1,:) 1920 final_bounds(2,:) = new_bounds(2,:) 1921 else 1922 final_bounds(1,:) = min( old_bounds(1,:), new_bounds(1,:) ) 1923 final_bounds(2,:) = max( old_bounds(2,:), new_bounds(2,:) ) 1924 end if 1925 1926 ! Find common section of old and new arrays 1927 common_bounds(1,:) = max( old_bounds(1,:), final_bounds(1,:) ) 1928 common_bounds(2,:) = min( old_bounds(2,:), final_bounds(2,:) ) 1929 end if 1930 1931 end if 1932 1933else 1934 ! Old array does not exist. Allocate new one 1935 NEEDS_ALLOC = .true. 1936 NEEDS_DEALLOC = .false. 1937 NEEDS_COPY = .false. 1938 final_bounds(1,:) = new_bounds(1,:) 1939 final_bounds(2,:) = new_bounds(2,:) 1940end if 1941 1942END SUBROUTINE options 1943 1944! ================================================================== 1945 1946SUBROUTINE alloc_count( delta_size, type, name, routine ) 1947 1948implicit none 1949 1950integer, intent(in) :: delta_size ! +/-size(array) 1951character, intent(in) :: type ! 'I' => integer 1952 ! 'E' => integer*8 1953 ! 'R' => real*4 1954 ! 'D' => real*8 1955 ! 'L' => logical 1956 ! 'S' => character (string) 1957character(len=*), optional, intent(in) :: name 1958character(len=*), optional, intent(in) :: routine 1959 1960character(len=32) :: aname, rname 1961character(len=1) :: memType, task 1962real(DP) :: delta_mem 1963logical :: newPeak 1964logical, save :: header_written = .false. 1965logical, save :: tree_nullified = .false. 1966integer :: memSize 1967 1968! Set routine name 1969if (present(routine)) then 1970 rname = routine 1971else 1972 rname = DEFAULT%routine 1973end if 1974 1975! Call siesta counting routine 'memory' 1976! Switched off in Aug.2009, and made 'memory' call alloc_count instead 1977! if (delta_size > 0) then 1978! task = 'A' 1979! else 1980! task = 'D' 1981! end if 1982! select case( type ) 1983! case ('R') ! Real --> single 1984! memType = 'S' 1985! memSize = abs(delta_size) 1986! case ('S') ! Character (string) --> integer/4 1987! memType = 'I' 1988! memSize = abs(delta_size) / type_mem('I') 1989! case default 1990! memType = type 1991! memSize = abs(delta_size) 1992! end select 1993! call memory( task, memType, memSize, trim(rname) ) 1994 1995 1996if (REPORT_LEVEL <= 0) return 1997 1998! Compound routine+array name 1999if (present(name) .and. present(routine)) then 2000 aname = trim(routine)//' '//name 2001else if (present(name) .and. DEFAULT%routine/=DEFAULT_ROUTINE) then 2002 aname = trim(DEFAULT%routine)//' '//name 2003else if (present(name)) then 2004 aname = name 2005else if (present(routine)) then 2006 aname = trim(routine)//' '//DEFAULT_NAME 2007else if (DEFAULT%routine/=DEFAULT_ROUTINE) then 2008 aname = trim(DEFAULT%routine)//' '//DEFAULT_NAME 2009else 2010 aname = DEFAULT_ROUTINE//' '//DEFAULT_NAME 2011end if 2012 2013MAX_LEN = max( MAX_LEN, len(trim(aname)) ) 2014 2015! Find memory increment and total allocated memory 2016delta_mem = delta_size * type_mem(type) 2017TOT_MEM = TOT_MEM + delta_mem 2018if (TOT_MEM > PEAK_MEM+0.5_dp) then 2019 newPeak = .true. 2020 PEAK_MEM = TOT_MEM 2021 PEAK_ARRAY = aname 2022 PEAK_ROUTINE = rname 2023! print'(/,a,f18.6),a,/)', 2024! 'alloc_count: Memory peak =', PEAK_MEM/MBYTE, ' Mbytes' 2025else 2026 newPeak = .false. 2027end if 2028 2029! Add/subtract/classify array memory 2030if (REPORT_LEVEL > 1) then 2031 if (.not.tree_nullified) then 2032 nullify(report_tree) 2033 tree_nullified = .true. 2034 end if 2035 call tree_add( report_tree, aname, delta_mem ) 2036 if (newPeak) call tree_peak( report_tree ) 2037end if 2038 2039! Print report, but only in node 0, as not all 2040! processors may follow the same route here 2041! The detail/extent of the report increses with the value of level: 2042! level=0 : no report at all (the default) 2043! level=1 : only total memory peak and where it occurred 2044! level=2 : detailed report created but printed only upon request 2045! level=3 : detailed report printed at every new memory peak 2046! level=4 : print every individual reallocation or deallocation 2047 2048 2049if (newPeak .and. (REPORT_LEVEL==1 .or. REPORT_LEVEL==3) .and. & 2050 node == 0) then 2051 call print_report(.false.) 2052end if 2053 2054if (REPORT_LEVEL == 4 .and. node == 0) then 2055 if (.not.header_written) then 2056 write(REPORT_UNIT,'(/,a7,9x,1x,a4,28x,1x,2a15)') & 2057 'Routine', 'Name', 'Incr. (MB)', 'Total (MB)' 2058 header_written = .true. 2059 end if 2060 write(REPORT_UNIT,'(a16,1x,a32,1x,2f15.6)') & 2061 rname, aname, delta_mem/MBYTE, TOT_MEM/MBYTE 2062end if 2063END SUBROUTINE alloc_count 2064 2065! ================================================================== 2066 2067INTEGER FUNCTION type_mem( var_type ) 2068! 2069! It is not clear that the sizes assumed are universal for 2070! non-Cray machines... 2071! 2072implicit none 2073character, intent(in) :: var_type 2074character(len=40) :: message 2075 2076select case( var_type ) 2077#ifdef OLD_CRAY 2078 case('I') 2079 type_mem = 8 2080 case('R') 2081 type_mem = 8 2082 case('L') 2083 type_mem = 8 2084#else 2085 case('I') 2086 type_mem = 4 2087 case('R') 2088 type_mem = 4 2089 case('L') 2090 type_mem = 4 2091#endif 2092case('E') 2093 type_mem = 8 2094case('D') 2095 type_mem = 8 2096case('S') 2097 type_mem = 1 2098case default 2099 write(message,"(2a)") & 2100 'alloc_count: ERROR: unknown type = ', var_type 2101 call die(trim(message)) 2102end select 2103 2104END FUNCTION type_mem 2105 2106! ================================================================== 2107 2108RECURSIVE SUBROUTINE tree_add( t, name, delta_mem ) 2109 2110implicit none 2111type(TREE), pointer :: t 2112character(len=*), intent(in) :: name 2113real(DP), intent(in) :: delta_mem 2114 2115logical, save :: warn_negative = .true. 2116 2117if (.not.associated(t)) then 2118 allocate( t ) 2119 t%name = name 2120 t%mem = delta_mem 2121 t%max = delta_mem 2122 t%peak = 0._dp 2123 nullify( t%left, t%right ) 2124else if (name == t%name) then 2125 t%mem = t%mem + delta_mem 2126 ! The abs is to handle the case of apparent de_allocs without re_allocs, 2127 ! caused by routine/name argument mismatches 2128 if (abs(t%mem) > abs(t%max)) t%max = t%mem 2129else if ( llt(name,t%name) ) then 2130 call tree_add( t%left, name, delta_mem ) 2131else 2132 call tree_add( t%right, name, delta_mem ) 2133end if 2134 2135if (warn_negative .and. t%mem<0._dp) then 2136 call parallel_init() ! Make sure that node and Nodes are initialized 2137 if (Node==0) then 2138 write(6,'(/,a,/,2a,/,a,f18.0,a)') & 2139 'WARNING: alloc-realloc-dealloc name mismatch', & 2140 ' Name: ', trim(name), & 2141 ' Size: ', t%mem, ' Bytes' 2142 if (Nodes>1) write(6,'(9x,a,i6)') 'Node:', Node 2143 write(6,'(9x,a)') 'Subsequent mismatches will not be reported' 2144 warn_negative = .false. ! Print this warning only once 2145 end if 2146end if 2147 2148END SUBROUTINE tree_add 2149 2150! ================================================================== 2151 2152RECURSIVE SUBROUTINE tree_peak( t ) 2153 2154implicit none 2155type(TREE), pointer :: t 2156 2157if (.not.associated(t)) return 2158 2159t%peak = t%mem 2160call tree_peak( t%left ) 2161call tree_peak( t%right ) 2162 2163END SUBROUTINE tree_peak 2164 2165! ================================================================== 2166 2167RECURSIVE SUBROUTINE tree_print( t ) 2168 2169implicit none 2170type(TREE), pointer :: t 2171 2172if (.not.associated(t)) return 2173 2174call tree_print( t%left ) 2175 2176if (abs(t%max) >= REPORT_THRESHOLD) then 2177 write(REPORT_UNIT,'(a,1x,3f15.6,f9.2)') & 2178 t%name(1:MAX_LEN), t%mem/MBYTE, t%max/MBYTE, t%peak/MBYTE, & 2179 100._dp * t%peak / (PEAK_MEM + tiny(PEAK_MEM) ) 2180end if 2181 2182call tree_print( t%right ) 2183 2184END SUBROUTINE tree_print 2185 2186! ================================================================== 2187 2188SUBROUTINE print_report(all) 2189 2190implicit none 2191 2192! Whether MPI reductions should be performed 2193! If, not then ensure that no MPI calls are performed 2194logical, intent(in), optional :: all 2195 2196character(len=80) :: string = 'Name' 2197character :: date*8, time*10, zone*5 2198integer :: iNode, peakNode 2199real(dp) :: maxPeak 2200real(dp),allocatable:: nodeMem(:), nodePeak(:) 2201logical :: lall 2202 2203#ifdef MPI 2204integer :: MPIerror 2205#endif 2206 2207! Enables parallel call of print-report 2208lall = .true. 2209if ( present(all) ) lall = all 2210 2211! Only if MPI should all be used 2212if ( lall ) then 2213! Make sure that variables node and Nodes are initialized 2214call parallel_init() 2215end if 2216 2217! Allocate and initialize two small arrays 2218allocate( nodeMem(0:Nodes-1), nodePeak(0:Nodes-1) ) 2219! initialize (in case all == .false.) 2220nodeMem = 0._dp 2221nodePeak = 0._dp 2222 2223! Initializations for Nodes=1 (serial case) 2224nodeMem(node) = TOT_MEM 2225nodePeak(node) = PEAK_MEM 2226peakNode = node 2227 2228! In parallel, find the memory values of all nodes 2229#ifdef MPI 2230if (Nodes > 1 .and. lall ) then 2231 ! Gather the present and peak memories of all nodes 2232 call MPI_AllGather( TOT_MEM, 1, MPI_double_precision, & 2233 nodeMem, 1, MPI_double_precision, & 2234 MPI_COMM_WORLD, MPIerror ) 2235 call MPI_AllGather( PEAK_MEM, 1, MPI_double_precision, & 2236 nodePeak, 1, MPI_double_precision, & 2237 MPI_COMM_WORLD, MPIerror ) 2238 ! Find the node with the highest peak of memory 2239 maxPeak = 0 2240 do iNode = 0,Nodes-1 2241 if (nodePeak(iNode) > maxPeak) then 2242 peakNode = iNode 2243 maxPeak = nodePeak(iNode) 2244 end if 2245 end do ! iNode 2246 ! Change the writing node for the peak-node information 2247 if (node==0 .and. peakNode/=0) close( REPORT_UNIT ) 2248 call MPI_Barrier( MPI_COMM_WORLD, MPIerror ) 2249 if (node==peakNode .and. peakNode/=0) then 2250 call io_assign(REPORT_UNIT) 2251 open( unit=REPORT_UNIT, file=REPORT_FILE, & 2252 status='unknown', position='append' ) 2253 end if 2254end if ! (Nodes>1) 2255#endif 2256 2257! The report is printed by the highest-peak node 2258if (node == peakNode) then 2259 2260 ! AG: Commented out to allow multiple batches of information 2261 ! if (REPORT_LEVEL < 4) rewind(REPORT_UNIT) 2262 2263 call date_and_time( date, time, zone ) 2264 2265 write(REPORT_UNIT,'(/,a,16a)') & 2266 'Allocation summary at ', & 2267 date(1:4),'/',date(5:6),'/',date(7:8),' ', & 2268 time(1:2),':',time(3:4),':',time(5:10),' ', & 2269 zone(1:3),':',zone(4:5) 2270 2271 if (Nodes > 1) then 2272 write(REPORT_UNIT,'(/,(a,f18.6,a))') & 2273 'Present memory all nodes : ', sum(nodeMem)/MBYTE, ' MB', & 2274 'Added peak mem all nodes : ', sum(nodePeak)/MBYTE, ' MB', & 2275 'Min peak memory in a node: ', minval(nodePeak)/MBYTE, ' MB', & 2276 'Max peak memory in a node: ', maxval(nodePeak)/MBYTE, ' MB' 2277! Impractical for many nodes: 2278! write(REPORT_UNIT,'(/,a,/,(i6,f12.6))') & 2279! 'Memory peaks of individual nodes (Mb):', & 2280! (iNode,nodePeak(iNode)/MBYTE,iNode=0,Nodes-1) 2281 write(REPORT_UNIT,'(/,a,i6)') & 2282 'Maximum peak of memory occurred in node:', peakNode 2283 end if 2284 2285 write(REPORT_UNIT,'(2(/,a,f18.6,a),/,2a,/,2a)') & 2286 'Present memory allocation: ', TOT_MEM/MBYTE, ' MB', & 2287 'Maximum memory allocation: ', PEAK_MEM/MBYTE, ' MB', & 2288 'Occurred after allocating: ', trim(PEAK_ARRAY), & 2289 'In routine: ', trim(PEAK_ROUTINE) 2290 2291 if (REPORT_LEVEL > 1) then 2292 if (REPORT_THRESHOLD > 0._dp) then 2293 write(REPORT_UNIT,'(/,a,f12.6,a,/,a,1x,3a15,a9)') & 2294 'Allocated sizes (in MByte) of arrays larger than ', & 2295 REPORT_THRESHOLD/MBYTE, ' MB:', & 2296 string(1:MAX_LEN), 'Present', 'Maximum', 'At peak', '%' 2297 else 2298 write(REPORT_UNIT,'(/,a,/,a,1x,3a15,a9)') & 2299 'Allocated array sizes (in MByte):', & 2300 string(1:MAX_LEN), 'Present', 'Maximum', 'At peak', '%' 2301 end if 2302 call tree_print( report_tree ) 2303 end if 2304 2305 ! Close file if not common IO node 2306 if ( node /= 0 ) call io_close(REPORT_UNIT) 2307end if ! (node == peakNode) 2308 2309! Change again the writing node for the rest of the report 2310#ifdef MPI 2311if (lall) call MPI_Barrier( MPI_COMM_WORLD, MPIerror ) 2312if (node==0 .and. peakNode/=0) & 2313 open( unit=REPORT_UNIT, file=REPORT_FILE, & 2314 status='unknown', position='append' ) 2315#endif 2316 2317deallocate( nodeMem, nodePeak ) 2318 2319END SUBROUTINE print_report 2320 2321! ================================================================== 2322 2323SUBROUTINE alloc_err( ierr, name, routine, bounds ) 2324 2325implicit none 2326 2327integer, intent(in) :: ierr 2328character(len=*), optional, intent(in) :: name 2329character(len=*), optional, intent(in) :: routine 2330integer, dimension(:,:), optional, intent(in) :: bounds 2331 2332integer i 2333 2334if (ierr/=0) then 2335 if (ionode) print*, 'alloc_err: allocate status error', ierr 2336 if (present(name).and.present(routine)) then 2337 if (ionode) print*, 'alloc_err: array ', name, & 2338 ' requested by ', routine 2339 elseif (present(name)) then 2340 if (ionode) print*, 'alloc_err: array ', name, & 2341 ' requested by unknown' 2342 elseif (present(routine)) then 2343 if (ionode) print* , 'alloc_err: array unknown', & 2344 ' requested by ', routine 2345 endif 2346 if (ionode.and.present(bounds)) & 2347 print '(a,i3,2i10)', ('alloc_err: dim, lbound, ubound:', & 2348 i,bounds(1,i),bounds(2,i), & 2349 i=1,size(bounds,dim=2)) 2350 2351 call die('alloc_err: allocate error') 2352end if 2353 2354END SUBROUTINE alloc_err 2355 2356! ================================================================== 2357 2358END MODULE alloc 2359