1 /* Implementation of the MINLOC intrinsic 2 Copyright (C) 2002-2019 Free Software Foundation, Inc. 3 Contributed by Paul Brook <paul@nowt.org> 4 5 This file is part of the GNU Fortran runtime library (libgfortran). 6 7 Libgfortran is free software; you can redistribute it and/or 8 modify it under the terms of the GNU General Public 9 License as published by the Free Software Foundation; either 10 version 3 of the License, or (at your option) any later version. 11 12 Libgfortran is distributed in the hope that it will be useful, 13 but WITHOUT ANY WARRANTY; without even the implied warranty of 14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 15 GNU General Public License for more details. 16 17 Under Section 7 of GPL version 3, you are granted additional 18 permissions described in the GCC Runtime Library Exception, version 19 3.1, as published by the Free Software Foundation. 20 21 You should have received a copy of the GNU General Public License and 22 a copy of the GCC Runtime Library Exception along with this program; 23 see the files COPYING3 and COPYING.RUNTIME respectively. If not, see 24 <http://www.gnu.org/licenses/>. */ 25 26 #include "libgfortran.h" 27 #include <assert.h> 28 29 30 #if defined (HAVE_GFC_REAL_4) && defined (HAVE_GFC_INTEGER_8) 31 32 #define HAVE_BACK_ARG 1 33 34 35 extern void minloc1_8_r4 (gfc_array_i8 * const restrict, 36 gfc_array_r4 * const restrict, const index_type * const restrict, GFC_LOGICAL_4 back); 37 export_proto(minloc1_8_r4); 38 39 void 40 minloc1_8_r4 (gfc_array_i8 * const restrict retarray, 41 gfc_array_r4 * const restrict array, 42 const index_type * const restrict pdim, GFC_LOGICAL_4 back) 43 { 44 index_type count[GFC_MAX_DIMENSIONS]; 45 index_type extent[GFC_MAX_DIMENSIONS]; 46 index_type sstride[GFC_MAX_DIMENSIONS]; 47 index_type dstride[GFC_MAX_DIMENSIONS]; 48 const GFC_REAL_4 * restrict base; 49 GFC_INTEGER_8 * restrict dest; 50 index_type rank; 51 index_type n; 52 index_type len; 53 index_type delta; 54 index_type dim; 55 int continue_loop; 56 57 /* Make dim zero based to avoid confusion. */ 58 rank = GFC_DESCRIPTOR_RANK (array) - 1; 59 dim = (*pdim) - 1; 60 61 if (unlikely (dim < 0 || dim > rank)) 62 { 63 runtime_error ("Dim argument incorrect in MINLOC intrinsic: " 64 "is %ld, should be between 1 and %ld", 65 (long int) dim + 1, (long int) rank + 1); 66 } 67 68 len = GFC_DESCRIPTOR_EXTENT(array,dim); 69 if (len < 0) 70 len = 0; 71 delta = GFC_DESCRIPTOR_STRIDE(array,dim); 72 73 for (n = 0; n < dim; n++) 74 { 75 sstride[n] = GFC_DESCRIPTOR_STRIDE(array,n); 76 extent[n] = GFC_DESCRIPTOR_EXTENT(array,n); 77 78 if (extent[n] < 0) 79 extent[n] = 0; 80 } 81 for (n = dim; n < rank; n++) 82 { 83 sstride[n] = GFC_DESCRIPTOR_STRIDE(array, n + 1); 84 extent[n] = GFC_DESCRIPTOR_EXTENT(array, n + 1); 85 86 if (extent[n] < 0) 87 extent[n] = 0; 88 } 89 90 if (retarray->base_addr == NULL) 91 { 92 size_t alloc_size, str; 93 94 for (n = 0; n < rank; n++) 95 { 96 if (n == 0) 97 str = 1; 98 else 99 str = GFC_DESCRIPTOR_STRIDE(retarray,n-1) * extent[n-1]; 100 101 GFC_DIMENSION_SET(retarray->dim[n], 0, extent[n] - 1, str); 102 103 } 104 105 retarray->offset = 0; 106 retarray->dtype.rank = rank; 107 108 alloc_size = GFC_DESCRIPTOR_STRIDE(retarray,rank-1) * extent[rank-1]; 109 110 retarray->base_addr = xmallocarray (alloc_size, sizeof (GFC_INTEGER_8)); 111 if (alloc_size == 0) 112 { 113 /* Make sure we have a zero-sized array. */ 114 GFC_DIMENSION_SET(retarray->dim[0], 0, -1, 1); 115 return; 116 117 } 118 } 119 else 120 { 121 if (rank != GFC_DESCRIPTOR_RANK (retarray)) 122 runtime_error ("rank of return array incorrect in" 123 " MINLOC intrinsic: is %ld, should be %ld", 124 (long int) (GFC_DESCRIPTOR_RANK (retarray)), 125 (long int) rank); 126 127 if (unlikely (compile_options.bounds_check)) 128 bounds_ifunction_return ((array_t *) retarray, extent, 129 "return value", "MINLOC"); 130 } 131 132 for (n = 0; n < rank; n++) 133 { 134 count[n] = 0; 135 dstride[n] = GFC_DESCRIPTOR_STRIDE(retarray,n); 136 if (extent[n] <= 0) 137 return; 138 } 139 140 base = array->base_addr; 141 dest = retarray->base_addr; 142 143 continue_loop = 1; 144 while (continue_loop) 145 { 146 const GFC_REAL_4 * restrict src; 147 GFC_INTEGER_8 result; 148 src = base; 149 { 150 151 GFC_REAL_4 minval; 152 #if defined (GFC_REAL_4_INFINITY) 153 minval = GFC_REAL_4_INFINITY; 154 #else 155 minval = GFC_REAL_4_HUGE; 156 #endif 157 result = 1; 158 if (len <= 0) 159 *dest = 0; 160 else 161 { 162 #if ! defined HAVE_BACK_ARG 163 for (n = 0; n < len; n++, src += delta) 164 { 165 #endif 166 167 #if defined (GFC_REAL_4_QUIET_NAN) 168 for (n = 0; n < len; n++, src += delta) 169 { 170 if (*src <= minval) 171 { 172 minval = *src; 173 result = (GFC_INTEGER_8)n + 1; 174 break; 175 } 176 } 177 #else 178 n = 0; 179 #endif 180 if (back) 181 for (; n < len; n++, src += delta) 182 { 183 if (unlikely (*src <= minval)) 184 { 185 minval = *src; 186 result = (GFC_INTEGER_8)n + 1; 187 } 188 } 189 else 190 for (; n < len; n++, src += delta) 191 { 192 if (unlikely (*src < minval)) 193 { 194 minval = *src; 195 result = (GFC_INTEGER_8) n + 1; 196 } 197 } 198 199 *dest = result; 200 } 201 } 202 /* Advance to the next element. */ 203 count[0]++; 204 base += sstride[0]; 205 dest += dstride[0]; 206 n = 0; 207 while (count[n] == extent[n]) 208 { 209 /* When we get to the end of a dimension, reset it and increment 210 the next dimension. */ 211 count[n] = 0; 212 /* We could precalculate these products, but this is a less 213 frequently used path so probably not worth it. */ 214 base -= sstride[n] * extent[n]; 215 dest -= dstride[n] * extent[n]; 216 n++; 217 if (n >= rank) 218 { 219 /* Break out of the loop. */ 220 continue_loop = 0; 221 break; 222 } 223 else 224 { 225 count[n]++; 226 base += sstride[n]; 227 dest += dstride[n]; 228 } 229 } 230 } 231 } 232 233 234 extern void mminloc1_8_r4 (gfc_array_i8 * const restrict, 235 gfc_array_r4 * const restrict, const index_type * const restrict, 236 gfc_array_l1 * const restrict, GFC_LOGICAL_4 back); 237 export_proto(mminloc1_8_r4); 238 239 void 240 mminloc1_8_r4 (gfc_array_i8 * const restrict retarray, 241 gfc_array_r4 * const restrict array, 242 const index_type * const restrict pdim, 243 gfc_array_l1 * const restrict mask, GFC_LOGICAL_4 back) 244 { 245 index_type count[GFC_MAX_DIMENSIONS]; 246 index_type extent[GFC_MAX_DIMENSIONS]; 247 index_type sstride[GFC_MAX_DIMENSIONS]; 248 index_type dstride[GFC_MAX_DIMENSIONS]; 249 index_type mstride[GFC_MAX_DIMENSIONS]; 250 GFC_INTEGER_8 * restrict dest; 251 const GFC_REAL_4 * restrict base; 252 const GFC_LOGICAL_1 * restrict mbase; 253 index_type rank; 254 index_type dim; 255 index_type n; 256 index_type len; 257 index_type delta; 258 index_type mdelta; 259 int mask_kind; 260 261 if (mask == NULL) 262 { 263 #ifdef HAVE_BACK_ARG 264 minloc1_8_r4 (retarray, array, pdim, back); 265 #else 266 minloc1_8_r4 (retarray, array, pdim); 267 #endif 268 return; 269 } 270 271 dim = (*pdim) - 1; 272 rank = GFC_DESCRIPTOR_RANK (array) - 1; 273 274 275 if (unlikely (dim < 0 || dim > rank)) 276 { 277 runtime_error ("Dim argument incorrect in MINLOC intrinsic: " 278 "is %ld, should be between 1 and %ld", 279 (long int) dim + 1, (long int) rank + 1); 280 } 281 282 len = GFC_DESCRIPTOR_EXTENT(array,dim); 283 if (len <= 0) 284 return; 285 286 mbase = mask->base_addr; 287 288 mask_kind = GFC_DESCRIPTOR_SIZE (mask); 289 290 if (mask_kind == 1 || mask_kind == 2 || mask_kind == 4 || mask_kind == 8 291 #ifdef HAVE_GFC_LOGICAL_16 292 || mask_kind == 16 293 #endif 294 ) 295 mbase = GFOR_POINTER_TO_L1 (mbase, mask_kind); 296 else 297 runtime_error ("Funny sized logical array"); 298 299 delta = GFC_DESCRIPTOR_STRIDE(array,dim); 300 mdelta = GFC_DESCRIPTOR_STRIDE_BYTES(mask,dim); 301 302 for (n = 0; n < dim; n++) 303 { 304 sstride[n] = GFC_DESCRIPTOR_STRIDE(array,n); 305 mstride[n] = GFC_DESCRIPTOR_STRIDE_BYTES(mask,n); 306 extent[n] = GFC_DESCRIPTOR_EXTENT(array,n); 307 308 if (extent[n] < 0) 309 extent[n] = 0; 310 311 } 312 for (n = dim; n < rank; n++) 313 { 314 sstride[n] = GFC_DESCRIPTOR_STRIDE(array,n + 1); 315 mstride[n] = GFC_DESCRIPTOR_STRIDE_BYTES(mask, n + 1); 316 extent[n] = GFC_DESCRIPTOR_EXTENT(array, n + 1); 317 318 if (extent[n] < 0) 319 extent[n] = 0; 320 } 321 322 if (retarray->base_addr == NULL) 323 { 324 size_t alloc_size, str; 325 326 for (n = 0; n < rank; n++) 327 { 328 if (n == 0) 329 str = 1; 330 else 331 str= GFC_DESCRIPTOR_STRIDE(retarray,n-1) * extent[n-1]; 332 333 GFC_DIMENSION_SET(retarray->dim[n], 0, extent[n] - 1, str); 334 335 } 336 337 alloc_size = GFC_DESCRIPTOR_STRIDE(retarray,rank-1) * extent[rank-1]; 338 339 retarray->offset = 0; 340 retarray->dtype.rank = rank; 341 342 if (alloc_size == 0) 343 { 344 /* Make sure we have a zero-sized array. */ 345 GFC_DIMENSION_SET(retarray->dim[0], 0, -1, 1); 346 return; 347 } 348 else 349 retarray->base_addr = xmallocarray (alloc_size, sizeof (GFC_INTEGER_8)); 350 351 } 352 else 353 { 354 if (rank != GFC_DESCRIPTOR_RANK (retarray)) 355 runtime_error ("rank of return array incorrect in MINLOC intrinsic"); 356 357 if (unlikely (compile_options.bounds_check)) 358 { 359 bounds_ifunction_return ((array_t *) retarray, extent, 360 "return value", "MINLOC"); 361 bounds_equal_extents ((array_t *) mask, (array_t *) array, 362 "MASK argument", "MINLOC"); 363 } 364 } 365 366 for (n = 0; n < rank; n++) 367 { 368 count[n] = 0; 369 dstride[n] = GFC_DESCRIPTOR_STRIDE(retarray,n); 370 if (extent[n] <= 0) 371 return; 372 } 373 374 dest = retarray->base_addr; 375 base = array->base_addr; 376 377 while (base) 378 { 379 const GFC_REAL_4 * restrict src; 380 const GFC_LOGICAL_1 * restrict msrc; 381 GFC_INTEGER_8 result; 382 src = base; 383 msrc = mbase; 384 { 385 386 GFC_REAL_4 minval; 387 #if defined (GFC_REAL_4_INFINITY) 388 minval = GFC_REAL_4_INFINITY; 389 #else 390 minval = GFC_REAL_4_HUGE; 391 #endif 392 #if defined (GFC_REAL_4_QUIET_NAN) 393 GFC_INTEGER_8 result2 = 0; 394 #endif 395 result = 0; 396 for (n = 0; n < len; n++, src += delta, msrc += mdelta) 397 { 398 399 if (*msrc) 400 { 401 #if defined (GFC_REAL_4_QUIET_NAN) 402 if (!result2) 403 result2 = (GFC_INTEGER_8)n + 1; 404 if (*src <= minval) 405 #endif 406 { 407 minval = *src; 408 result = (GFC_INTEGER_8)n + 1; 409 break; 410 } 411 } 412 } 413 #if defined (GFC_REAL_4_QUIET_NAN) 414 if (unlikely (n >= len)) 415 result = result2; 416 else 417 #endif 418 if (back) 419 for (; n < len; n++, src += delta, msrc += mdelta) 420 { 421 if (*msrc && unlikely (*src <= minval)) 422 { 423 minval = *src; 424 result = (GFC_INTEGER_8)n + 1; 425 } 426 } 427 else 428 for (; n < len; n++, src += delta, msrc += mdelta) 429 { 430 if (*msrc && unlikely (*src < minval)) 431 { 432 minval = *src; 433 result = (GFC_INTEGER_8) n + 1; 434 } 435 } 436 *dest = result; 437 } 438 /* Advance to the next element. */ 439 count[0]++; 440 base += sstride[0]; 441 mbase += mstride[0]; 442 dest += dstride[0]; 443 n = 0; 444 while (count[n] == extent[n]) 445 { 446 /* When we get to the end of a dimension, reset it and increment 447 the next dimension. */ 448 count[n] = 0; 449 /* We could precalculate these products, but this is a less 450 frequently used path so probably not worth it. */ 451 base -= sstride[n] * extent[n]; 452 mbase -= mstride[n] * extent[n]; 453 dest -= dstride[n] * extent[n]; 454 n++; 455 if (n >= rank) 456 { 457 /* Break out of the loop. */ 458 base = NULL; 459 break; 460 } 461 else 462 { 463 count[n]++; 464 base += sstride[n]; 465 mbase += mstride[n]; 466 dest += dstride[n]; 467 } 468 } 469 } 470 } 471 472 473 extern void sminloc1_8_r4 (gfc_array_i8 * const restrict, 474 gfc_array_r4 * const restrict, const index_type * const restrict, 475 GFC_LOGICAL_4 *, GFC_LOGICAL_4 back); 476 export_proto(sminloc1_8_r4); 477 478 void 479 sminloc1_8_r4 (gfc_array_i8 * const restrict retarray, 480 gfc_array_r4 * const restrict array, 481 const index_type * const restrict pdim, 482 GFC_LOGICAL_4 * mask, GFC_LOGICAL_4 back) 483 { 484 index_type count[GFC_MAX_DIMENSIONS]; 485 index_type extent[GFC_MAX_DIMENSIONS]; 486 index_type dstride[GFC_MAX_DIMENSIONS]; 487 GFC_INTEGER_8 * restrict dest; 488 index_type rank; 489 index_type n; 490 index_type dim; 491 492 493 if (mask == NULL || *mask) 494 { 495 #ifdef HAVE_BACK_ARG 496 minloc1_8_r4 (retarray, array, pdim, back); 497 #else 498 minloc1_8_r4 (retarray, array, pdim); 499 #endif 500 return; 501 } 502 /* Make dim zero based to avoid confusion. */ 503 dim = (*pdim) - 1; 504 rank = GFC_DESCRIPTOR_RANK (array) - 1; 505 506 if (unlikely (dim < 0 || dim > rank)) 507 { 508 runtime_error ("Dim argument incorrect in MINLOC intrinsic: " 509 "is %ld, should be between 1 and %ld", 510 (long int) dim + 1, (long int) rank + 1); 511 } 512 513 for (n = 0; n < dim; n++) 514 { 515 extent[n] = GFC_DESCRIPTOR_EXTENT(array,n); 516 517 if (extent[n] <= 0) 518 extent[n] = 0; 519 } 520 521 for (n = dim; n < rank; n++) 522 { 523 extent[n] = 524 GFC_DESCRIPTOR_EXTENT(array,n + 1); 525 526 if (extent[n] <= 0) 527 extent[n] = 0; 528 } 529 530 if (retarray->base_addr == NULL) 531 { 532 size_t alloc_size, str; 533 534 for (n = 0; n < rank; n++) 535 { 536 if (n == 0) 537 str = 1; 538 else 539 str = GFC_DESCRIPTOR_STRIDE(retarray,n-1) * extent[n-1]; 540 541 GFC_DIMENSION_SET(retarray->dim[n], 0, extent[n] - 1, str); 542 543 } 544 545 retarray->offset = 0; 546 retarray->dtype.rank = rank; 547 548 alloc_size = GFC_DESCRIPTOR_STRIDE(retarray,rank-1) * extent[rank-1]; 549 550 if (alloc_size == 0) 551 { 552 /* Make sure we have a zero-sized array. */ 553 GFC_DIMENSION_SET(retarray->dim[0], 0, -1, 1); 554 return; 555 } 556 else 557 retarray->base_addr = xmallocarray (alloc_size, sizeof (GFC_INTEGER_8)); 558 } 559 else 560 { 561 if (rank != GFC_DESCRIPTOR_RANK (retarray)) 562 runtime_error ("rank of return array incorrect in" 563 " MINLOC intrinsic: is %ld, should be %ld", 564 (long int) (GFC_DESCRIPTOR_RANK (retarray)), 565 (long int) rank); 566 567 if (unlikely (compile_options.bounds_check)) 568 { 569 for (n=0; n < rank; n++) 570 { 571 index_type ret_extent; 572 573 ret_extent = GFC_DESCRIPTOR_EXTENT(retarray,n); 574 if (extent[n] != ret_extent) 575 runtime_error ("Incorrect extent in return value of" 576 " MINLOC intrinsic in dimension %ld:" 577 " is %ld, should be %ld", (long int) n + 1, 578 (long int) ret_extent, (long int) extent[n]); 579 } 580 } 581 } 582 583 for (n = 0; n < rank; n++) 584 { 585 count[n] = 0; 586 dstride[n] = GFC_DESCRIPTOR_STRIDE(retarray,n); 587 } 588 589 dest = retarray->base_addr; 590 591 while(1) 592 { 593 *dest = 0; 594 count[0]++; 595 dest += dstride[0]; 596 n = 0; 597 while (count[n] == extent[n]) 598 { 599 /* When we get to the end of a dimension, reset it and increment 600 the next dimension. */ 601 count[n] = 0; 602 /* We could precalculate these products, but this is a less 603 frequently used path so probably not worth it. */ 604 dest -= dstride[n] * extent[n]; 605 n++; 606 if (n >= rank) 607 return; 608 else 609 { 610 count[n]++; 611 dest += dstride[n]; 612 } 613 } 614 } 615 } 616 617 #endif 618