minloc1_8_r16.c revision 1.1.1.3
1/* Implementation of the MINLOC intrinsic 2 Copyright (C) 2002-2022 Free Software Foundation, Inc. 3 Contributed by Paul Brook <paul@nowt.org> 4 5This file is part of the GNU Fortran runtime library (libgfortran). 6 7Libgfortran is free software; you can redistribute it and/or 8modify it under the terms of the GNU General Public 9License as published by the Free Software Foundation; either 10version 3 of the License, or (at your option) any later version. 11 12Libgfortran is distributed in the hope that it will be useful, 13but WITHOUT ANY WARRANTY; without even the implied warranty of 14MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 15GNU General Public License for more details. 16 17Under Section 7 of GPL version 3, you are granted additional 18permissions described in the GCC Runtime Library Exception, version 193.1, as published by the Free Software Foundation. 20 21You should have received a copy of the GNU General Public License and 22a copy of the GCC Runtime Library Exception along with this program; 23see 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_16) && defined (HAVE_GFC_INTEGER_8) 31 32#define HAVE_BACK_ARG 1 33 34 35extern void minloc1_8_r16 (gfc_array_i8 * const restrict, 36 gfc_array_r16 * const restrict, const index_type * const restrict, GFC_LOGICAL_4 back); 37export_proto(minloc1_8_r16); 38 39void 40minloc1_8_r16 (gfc_array_i8 * const restrict retarray, 41 gfc_array_r16 * 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_16 * 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_16 * restrict src; 147 GFC_INTEGER_8 result; 148 src = base; 149 { 150 151 GFC_REAL_16 minval; 152#if defined (GFC_REAL_16_INFINITY) 153 minval = GFC_REAL_16_INFINITY; 154#else 155 minval = GFC_REAL_16_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_16_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 234extern void mminloc1_8_r16 (gfc_array_i8 * const restrict, 235 gfc_array_r16 * const restrict, const index_type * const restrict, 236 gfc_array_l1 * const restrict, GFC_LOGICAL_4 back); 237export_proto(mminloc1_8_r16); 238 239void 240mminloc1_8_r16 (gfc_array_i8 * const restrict retarray, 241 gfc_array_r16 * 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_16 * 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_r16 (retarray, array, pdim, back); 265#else 266 minloc1_8_r16 (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_16 * 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_16 minval; 387#if defined (GFC_REAL_16_INFINITY) 388 minval = GFC_REAL_16_INFINITY; 389#else 390 minval = GFC_REAL_16_HUGE; 391#endif 392#if defined (GFC_REAL_16_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_16_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_16_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 473extern void sminloc1_8_r16 (gfc_array_i8 * const restrict, 474 gfc_array_r16 * const restrict, const index_type * const restrict, 475 GFC_LOGICAL_4 *, GFC_LOGICAL_4 back); 476export_proto(sminloc1_8_r16); 477 478void 479sminloc1_8_r16 (gfc_array_i8 * const restrict retarray, 480 gfc_array_r16 * 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_r16 (retarray, array, pdim, back); 497#else 498 minloc1_8_r16 (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