197403Sobrien/* Implementation of the PRODUCT intrinsic 297403Sobrien Copyright (C) 2002-2020 Free Software Foundation, Inc. 397403Sobrien Contributed by Paul Brook <paul@nowt.org> 497403Sobrien 597403SobrienThis file is part of the GNU Fortran runtime library (libgfortran). 697403Sobrien 797403SobrienLibgfortran is free software; you can redistribute it and/or 897403Sobrienmodify it under the terms of the GNU General Public 997403SobrienLicense as published by the Free Software Foundation; either 1097403Sobrienversion 3 of the License, or (at your option) any later version. 1197403Sobrien 1297403SobrienLibgfortran is distributed in the hope that it will be useful, 1397403Sobrienbut WITHOUT ANY WARRANTY; without even the implied warranty of 1497403SobrienMERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 1597403SobrienGNU General Public License for more details. 1697403Sobrien 1797403SobrienUnder Section 7 of GPL version 3, you are granted additional 1897403Sobrienpermissions described in the GCC Runtime Library Exception, version 1997403Sobrien3.1, as published by the Free Software Foundation. 2097403Sobrien 2197403SobrienYou should have received a copy of the GNU General Public License and 2297403Sobriena copy of the GCC Runtime Library Exception along with this program; 2397403Sobriensee the files COPYING3 and COPYING.RUNTIME respectively. If not, see 2497403Sobrien<http://www.gnu.org/licenses/>. */ 2597403Sobrien 2697403Sobrien#include "libgfortran.h" 2797403Sobrien 2897403Sobrien 2997403Sobrien#if defined (HAVE_GFC_INTEGER_4) && defined (HAVE_GFC_INTEGER_4) 3097403Sobrien 3197403Sobrien 3297403Sobrienextern void product_i4 (gfc_array_i4 * const restrict, 3397403Sobrien gfc_array_i4 * const restrict, const index_type * const restrict); 3497403Sobrienexport_proto(product_i4); 3597403Sobrien 3697403Sobrienvoid 3797403Sobrienproduct_i4 (gfc_array_i4 * const restrict retarray, 3897403Sobrien gfc_array_i4 * const restrict array, 3997403Sobrien const index_type * const restrict pdim) 4097403Sobrien{ 4197403Sobrien index_type count[GFC_MAX_DIMENSIONS]; 4297403Sobrien index_type extent[GFC_MAX_DIMENSIONS]; 4397403Sobrien index_type sstride[GFC_MAX_DIMENSIONS]; 4497403Sobrien index_type dstride[GFC_MAX_DIMENSIONS]; 4597403Sobrien const GFC_INTEGER_4 * restrict base; 4697403Sobrien GFC_INTEGER_4 * restrict dest; 4797403Sobrien index_type rank; 4897403Sobrien index_type n; 4997403Sobrien index_type len; 5097403Sobrien index_type delta; 5197403Sobrien index_type dim; 5297403Sobrien int continue_loop; 5397403Sobrien 5497403Sobrien /* Make dim zero based to avoid confusion. */ 5597403Sobrien rank = GFC_DESCRIPTOR_RANK (array) - 1; 5697403Sobrien dim = (*pdim) - 1; 5797403Sobrien 5897403Sobrien if (unlikely (dim < 0 || dim > rank)) 5997403Sobrien { 6097403Sobrien runtime_error ("Dim argument incorrect in PRODUCT intrinsic: " 6197403Sobrien "is %ld, should be between 1 and %ld", 6297403Sobrien (long int) dim + 1, (long int) rank + 1); 6397403Sobrien } 6497403Sobrien 6597403Sobrien len = GFC_DESCRIPTOR_EXTENT(array,dim); 6697403Sobrien if (len < 0) 6797403Sobrien len = 0; 6897403Sobrien delta = GFC_DESCRIPTOR_STRIDE(array,dim); 6997403Sobrien 7097403Sobrien for (n = 0; n < dim; n++) 7197403Sobrien { 7297403Sobrien sstride[n] = GFC_DESCRIPTOR_STRIDE(array,n); 7397403Sobrien extent[n] = GFC_DESCRIPTOR_EXTENT(array,n); 7497403Sobrien 7597403Sobrien if (extent[n] < 0) 7697403Sobrien extent[n] = 0; 7797403Sobrien } 7897403Sobrien for (n = dim; n < rank; n++) 79110614Skan { 8097403Sobrien sstride[n] = GFC_DESCRIPTOR_STRIDE(array, n + 1); 8197403Sobrien extent[n] = GFC_DESCRIPTOR_EXTENT(array, n + 1); 8297403Sobrien 8397403Sobrien if (extent[n] < 0) 8497403Sobrien extent[n] = 0; 8597403Sobrien } 8697403Sobrien 8797403Sobrien if (retarray->base_addr == NULL) 8897403Sobrien { 8997403Sobrien size_t alloc_size, str; 9097403Sobrien 9197403Sobrien for (n = 0; n < rank; n++) 9297403Sobrien { 9397403Sobrien if (n == 0) 9497403Sobrien str = 1; 9597403Sobrien else 9697403Sobrien str = GFC_DESCRIPTOR_STRIDE(retarray,n-1) * extent[n-1]; 9797403Sobrien 9897403Sobrien GFC_DIMENSION_SET(retarray->dim[n], 0, extent[n] - 1, str); 9997403Sobrien 10097403Sobrien } 10197403Sobrien 10297403Sobrien retarray->offset = 0; 10397403Sobrien retarray->dtype.rank = rank; 10497403Sobrien 10597403Sobrien alloc_size = GFC_DESCRIPTOR_STRIDE(retarray,rank-1) * extent[rank-1]; 10697403Sobrien 10797403Sobrien retarray->base_addr = xmallocarray (alloc_size, sizeof (GFC_INTEGER_4)); 10897403Sobrien if (alloc_size == 0) 10997403Sobrien { 11097403Sobrien /* Make sure we have a zero-sized array. */ 11197403Sobrien GFC_DIMENSION_SET(retarray->dim[0], 0, -1, 1); 11297403Sobrien return; 11397403Sobrien 114102782Skan } 11597403Sobrien } 11697403Sobrien else 11797403Sobrien { 11897403Sobrien if (rank != GFC_DESCRIPTOR_RANK (retarray)) 11997403Sobrien runtime_error ("rank of return array incorrect in" 12097403Sobrien " PRODUCT intrinsic: is %ld, should be %ld", 12197403Sobrien (long int) (GFC_DESCRIPTOR_RANK (retarray)), 12297403Sobrien (long int) rank); 12397403Sobrien 12497403Sobrien if (unlikely (compile_options.bounds_check)) 12597403Sobrien bounds_ifunction_return ((array_t *) retarray, extent, 12697403Sobrien "return value", "PRODUCT"); 12797403Sobrien } 12897403Sobrien 12997403Sobrien for (n = 0; n < rank; n++) 13097403Sobrien { 13197403Sobrien count[n] = 0; 13297403Sobrien dstride[n] = GFC_DESCRIPTOR_STRIDE(retarray,n); 13397403Sobrien if (extent[n] <= 0) 13497403Sobrien return; 13597403Sobrien } 13697403Sobrien 13797403Sobrien base = array->base_addr; 13897403Sobrien dest = retarray->base_addr; 13997403Sobrien 14097403Sobrien continue_loop = 1; 14197403Sobrien while (continue_loop) 14297403Sobrien { 14397403Sobrien const GFC_INTEGER_4 * restrict src; 14497403Sobrien GFC_INTEGER_4 result; 14597403Sobrien src = base; 14697403Sobrien { 14797403Sobrien 14897403Sobrien result = 1; 14997403Sobrien if (len <= 0) 15097403Sobrien *dest = 1; 15197403Sobrien else 15297403Sobrien { 153103447Skan#if ! defined HAVE_BACK_ARG 154103447Skan for (n = 0; n < len; n++, src += delta) 155103447Skan { 15697403Sobrien#endif 15797403Sobrien 15897403Sobrien result *= *src; 15997403Sobrien } 16097403Sobrien 16197403Sobrien *dest = result; 16297403Sobrien } 16397403Sobrien } 16497403Sobrien /* Advance to the next element. */ 16597403Sobrien count[0]++; 16697403Sobrien base += sstride[0]; 16797403Sobrien dest += dstride[0]; 16897403Sobrien n = 0; 16997403Sobrien while (count[n] == extent[n]) 17097403Sobrien { 171103447Skan /* When we get to the end of a dimension, reset it and increment 17297403Sobrien the next dimension. */ 17397403Sobrien count[n] = 0; 17497403Sobrien /* We could precalculate these products, but this is a less 17597403Sobrien frequently used path so probably not worth it. */ 17697403Sobrien base -= sstride[n] * extent[n]; 17797403Sobrien dest -= dstride[n] * extent[n]; 17897403Sobrien n++; 17997403Sobrien if (n >= rank) 18097403Sobrien { 18197403Sobrien /* Break out of the loop. */ 18297403Sobrien continue_loop = 0; 18397403Sobrien break; 18497403Sobrien } 18597403Sobrien else 18697403Sobrien { 18797403Sobrien count[n]++; 18897403Sobrien base += sstride[n]; 18997403Sobrien dest += dstride[n]; 19097403Sobrien } 19197403Sobrien } 19297403Sobrien } 19397403Sobrien} 19497403Sobrien 19597403Sobrien 19697403Sobrienextern void mproduct_i4 (gfc_array_i4 * const restrict, 19797403Sobrien gfc_array_i4 * const restrict, const index_type * const restrict, 19897403Sobrien gfc_array_l1 * const restrict); 19997403Sobrienexport_proto(mproduct_i4); 20097403Sobrien 20197403Sobrienvoid 20297403Sobrienmproduct_i4 (gfc_array_i4 * const restrict retarray, 20397403Sobrien gfc_array_i4 * const restrict array, 20497403Sobrien const index_type * const restrict pdim, 20597403Sobrien gfc_array_l1 * const restrict mask) 20697403Sobrien{ 20797403Sobrien index_type count[GFC_MAX_DIMENSIONS]; 20897403Sobrien index_type extent[GFC_MAX_DIMENSIONS]; 20997403Sobrien index_type sstride[GFC_MAX_DIMENSIONS]; 21097403Sobrien index_type dstride[GFC_MAX_DIMENSIONS]; 21197403Sobrien index_type mstride[GFC_MAX_DIMENSIONS]; 21297403Sobrien GFC_INTEGER_4 * restrict dest; 21397403Sobrien const GFC_INTEGER_4 * restrict base; 21497403Sobrien const GFC_LOGICAL_1 * restrict mbase; 21597403Sobrien index_type rank; 21697403Sobrien index_type dim; 21797403Sobrien index_type n; 21897403Sobrien index_type len; 21997403Sobrien index_type delta; 22097403Sobrien index_type mdelta; 22197403Sobrien int mask_kind; 22297403Sobrien 22397403Sobrien if (mask == NULL) 22497403Sobrien { 22597403Sobrien#ifdef HAVE_BACK_ARG 22697403Sobrien product_i4 (retarray, array, pdim, back); 22797403Sobrien#else 22897403Sobrien product_i4 (retarray, array, pdim); 22997403Sobrien#endif 23097403Sobrien return; 23197403Sobrien } 23297403Sobrien 23397403Sobrien dim = (*pdim) - 1; 23497403Sobrien rank = GFC_DESCRIPTOR_RANK (array) - 1; 23597403Sobrien 23697403Sobrien 23797403Sobrien if (unlikely (dim < 0 || dim > rank)) 23897403Sobrien { 23997403Sobrien runtime_error ("Dim argument incorrect in PRODUCT intrinsic: " 24097403Sobrien "is %ld, should be between 1 and %ld", 24197403Sobrien (long int) dim + 1, (long int) rank + 1); 24297403Sobrien } 24397403Sobrien 24497403Sobrien len = GFC_DESCRIPTOR_EXTENT(array,dim); 24597403Sobrien if (len <= 0) 24697403Sobrien return; 24797403Sobrien 24897403Sobrien mbase = mask->base_addr; 24997403Sobrien 25097403Sobrien mask_kind = GFC_DESCRIPTOR_SIZE (mask); 25197403Sobrien 25297403Sobrien if (mask_kind == 1 || mask_kind == 2 || mask_kind == 4 || mask_kind == 8 25397403Sobrien#ifdef HAVE_GFC_LOGICAL_16 25497403Sobrien || mask_kind == 16 25597403Sobrien#endif 25697403Sobrien ) 25797403Sobrien mbase = GFOR_POINTER_TO_L1 (mbase, mask_kind); 25897403Sobrien else 25997403Sobrien runtime_error ("Funny sized logical array"); 26097403Sobrien 26197403Sobrien delta = GFC_DESCRIPTOR_STRIDE(array,dim); 26297403Sobrien mdelta = GFC_DESCRIPTOR_STRIDE_BYTES(mask,dim); 26397403Sobrien 26497403Sobrien for (n = 0; n < dim; n++) 26597403Sobrien { 26697403Sobrien sstride[n] = GFC_DESCRIPTOR_STRIDE(array,n); 26797403Sobrien mstride[n] = GFC_DESCRIPTOR_STRIDE_BYTES(mask,n); 26897403Sobrien extent[n] = GFC_DESCRIPTOR_EXTENT(array,n); 26997403Sobrien 27097403Sobrien if (extent[n] < 0) 27197403Sobrien extent[n] = 0; 27297403Sobrien 27397403Sobrien } 27497403Sobrien for (n = dim; n < rank; n++) 27597403Sobrien { 27697403Sobrien sstride[n] = GFC_DESCRIPTOR_STRIDE(array,n + 1); 27797403Sobrien mstride[n] = GFC_DESCRIPTOR_STRIDE_BYTES(mask, n + 1); 27897403Sobrien extent[n] = GFC_DESCRIPTOR_EXTENT(array, n + 1); 27997403Sobrien 28097403Sobrien if (extent[n] < 0) 28197403Sobrien extent[n] = 0; 28297403Sobrien } 28397403Sobrien 28497403Sobrien if (retarray->base_addr == NULL) 28597403Sobrien { 28697403Sobrien size_t alloc_size, str; 28797403Sobrien 28897403Sobrien for (n = 0; n < rank; n++) 28997403Sobrien { 29097403Sobrien if (n == 0) 29197403Sobrien str = 1; 29297403Sobrien else 29397403Sobrien str= GFC_DESCRIPTOR_STRIDE(retarray,n-1) * extent[n-1]; 29497403Sobrien 29597403Sobrien GFC_DIMENSION_SET(retarray->dim[n], 0, extent[n] - 1, str); 29697403Sobrien 29797403Sobrien } 29897403Sobrien 29997403Sobrien alloc_size = GFC_DESCRIPTOR_STRIDE(retarray,rank-1) * extent[rank-1]; 30097403Sobrien 30197403Sobrien retarray->offset = 0; 30297403Sobrien retarray->dtype.rank = rank; 30397403Sobrien 30497403Sobrien if (alloc_size == 0) 30597403Sobrien { 30697403Sobrien /* Make sure we have a zero-sized array. */ 30797403Sobrien GFC_DIMENSION_SET(retarray->dim[0], 0, -1, 1); 30897403Sobrien return; 30997403Sobrien } 31097403Sobrien else 31197403Sobrien retarray->base_addr = xmallocarray (alloc_size, sizeof (GFC_INTEGER_4)); 31297403Sobrien 31397403Sobrien } 31497403Sobrien else 31597403Sobrien { 31697403Sobrien if (rank != GFC_DESCRIPTOR_RANK (retarray)) 31797403Sobrien runtime_error ("rank of return array incorrect in PRODUCT intrinsic"); 31897403Sobrien 31997403Sobrien if (unlikely (compile_options.bounds_check)) 32097403Sobrien { 32197403Sobrien bounds_ifunction_return ((array_t *) retarray, extent, 32297403Sobrien "return value", "PRODUCT"); 32397403Sobrien bounds_equal_extents ((array_t *) mask, (array_t *) array, 32497403Sobrien "MASK argument", "PRODUCT"); 32597403Sobrien } 32697403Sobrien } 32797403Sobrien 32897403Sobrien for (n = 0; n < rank; n++) 32997403Sobrien { 33097403Sobrien count[n] = 0; 33197403Sobrien dstride[n] = GFC_DESCRIPTOR_STRIDE(retarray,n); 33297403Sobrien if (extent[n] <= 0) 33397403Sobrien return; 33497403Sobrien } 33597403Sobrien 33697403Sobrien dest = retarray->base_addr; 33797403Sobrien base = array->base_addr; 33897403Sobrien 33997403Sobrien while (base) 34097403Sobrien { 34197403Sobrien const GFC_INTEGER_4 * restrict src; 34297403Sobrien const GFC_LOGICAL_1 * restrict msrc; 34397403Sobrien GFC_INTEGER_4 result; 34497403Sobrien src = base; 34597403Sobrien msrc = mbase; 34697403Sobrien { 34797403Sobrien 34897403Sobrien result = 1; 34997403Sobrien for (n = 0; n < len; n++, src += delta, msrc += mdelta) 35097403Sobrien { 35197403Sobrien 35297403Sobrien if (*msrc) 35397403Sobrien result *= *src; 35497403Sobrien } 35597403Sobrien *dest = result; 35697403Sobrien } 35797403Sobrien /* Advance to the next element. */ 35897403Sobrien count[0]++; 35997403Sobrien base += sstride[0]; 36097403Sobrien mbase += mstride[0]; 36197403Sobrien dest += dstride[0]; 36297403Sobrien n = 0; 36397403Sobrien while (count[n] == extent[n]) 36497403Sobrien { 36597403Sobrien /* When we get to the end of a dimension, reset it and increment 36697403Sobrien the next dimension. */ 36797403Sobrien count[n] = 0; 36897403Sobrien /* We could precalculate these products, but this is a less 36997403Sobrien frequently used path so probably not worth it. */ 37097403Sobrien base -= sstride[n] * extent[n]; 37197403Sobrien mbase -= mstride[n] * extent[n]; 37297403Sobrien dest -= dstride[n] * extent[n]; 373 n++; 374 if (n >= rank) 375 { 376 /* Break out of the loop. */ 377 base = NULL; 378 break; 379 } 380 else 381 { 382 count[n]++; 383 base += sstride[n]; 384 mbase += mstride[n]; 385 dest += dstride[n]; 386 } 387 } 388 } 389} 390 391 392extern void sproduct_i4 (gfc_array_i4 * const restrict, 393 gfc_array_i4 * const restrict, const index_type * const restrict, 394 GFC_LOGICAL_4 *); 395export_proto(sproduct_i4); 396 397void 398sproduct_i4 (gfc_array_i4 * const restrict retarray, 399 gfc_array_i4 * const restrict array, 400 const index_type * const restrict pdim, 401 GFC_LOGICAL_4 * mask) 402{ 403 index_type count[GFC_MAX_DIMENSIONS]; 404 index_type extent[GFC_MAX_DIMENSIONS]; 405 index_type dstride[GFC_MAX_DIMENSIONS]; 406 GFC_INTEGER_4 * restrict dest; 407 index_type rank; 408 index_type n; 409 index_type dim; 410 411 412 if (mask == NULL || *mask) 413 { 414#ifdef HAVE_BACK_ARG 415 product_i4 (retarray, array, pdim, back); 416#else 417 product_i4 (retarray, array, pdim); 418#endif 419 return; 420 } 421 /* Make dim zero based to avoid confusion. */ 422 dim = (*pdim) - 1; 423 rank = GFC_DESCRIPTOR_RANK (array) - 1; 424 425 if (unlikely (dim < 0 || dim > rank)) 426 { 427 runtime_error ("Dim argument incorrect in PRODUCT intrinsic: " 428 "is %ld, should be between 1 and %ld", 429 (long int) dim + 1, (long int) rank + 1); 430 } 431 432 for (n = 0; n < dim; n++) 433 { 434 extent[n] = GFC_DESCRIPTOR_EXTENT(array,n); 435 436 if (extent[n] <= 0) 437 extent[n] = 0; 438 } 439 440 for (n = dim; n < rank; n++) 441 { 442 extent[n] = 443 GFC_DESCRIPTOR_EXTENT(array,n + 1); 444 445 if (extent[n] <= 0) 446 extent[n] = 0; 447 } 448 449 if (retarray->base_addr == NULL) 450 { 451 size_t alloc_size, str; 452 453 for (n = 0; n < rank; n++) 454 { 455 if (n == 0) 456 str = 1; 457 else 458 str = GFC_DESCRIPTOR_STRIDE(retarray,n-1) * extent[n-1]; 459 460 GFC_DIMENSION_SET(retarray->dim[n], 0, extent[n] - 1, str); 461 462 } 463 464 retarray->offset = 0; 465 retarray->dtype.rank = rank; 466 467 alloc_size = GFC_DESCRIPTOR_STRIDE(retarray,rank-1) * extent[rank-1]; 468 469 if (alloc_size == 0) 470 { 471 /* Make sure we have a zero-sized array. */ 472 GFC_DIMENSION_SET(retarray->dim[0], 0, -1, 1); 473 return; 474 } 475 else 476 retarray->base_addr = xmallocarray (alloc_size, sizeof (GFC_INTEGER_4)); 477 } 478 else 479 { 480 if (rank != GFC_DESCRIPTOR_RANK (retarray)) 481 runtime_error ("rank of return array incorrect in" 482 " PRODUCT intrinsic: is %ld, should be %ld", 483 (long int) (GFC_DESCRIPTOR_RANK (retarray)), 484 (long int) rank); 485 486 if (unlikely (compile_options.bounds_check)) 487 { 488 for (n=0; n < rank; n++) 489 { 490 index_type ret_extent; 491 492 ret_extent = GFC_DESCRIPTOR_EXTENT(retarray,n); 493 if (extent[n] != ret_extent) 494 runtime_error ("Incorrect extent in return value of" 495 " PRODUCT intrinsic in dimension %ld:" 496 " is %ld, should be %ld", (long int) n + 1, 497 (long int) ret_extent, (long int) extent[n]); 498 } 499 } 500 } 501 502 for (n = 0; n < rank; n++) 503 { 504 count[n] = 0; 505 dstride[n] = GFC_DESCRIPTOR_STRIDE(retarray,n); 506 } 507 508 dest = retarray->base_addr; 509 510 while(1) 511 { 512 *dest = 1; 513 count[0]++; 514 dest += dstride[0]; 515 n = 0; 516 while (count[n] == extent[n]) 517 { 518 /* When we get to the end of a dimension, reset it and increment 519 the next dimension. */ 520 count[n] = 0; 521 /* We could precalculate these products, but this is a less 522 frequently used path so probably not worth it. */ 523 dest -= dstride[n] * extent[n]; 524 n++; 525 if (n >= rank) 526 return; 527 else 528 { 529 count[n]++; 530 dest += dstride[n]; 531 } 532 } 533 } 534} 535 536#endif 537