minval_i2.c 12.8 KB
Newer Older
1
/* Implementation of the MINVAL intrinsic
2
   Copyright (C) 2002-2013 Free Software Foundation, Inc.
3 4
   Contributed by Paul Brook <paul@nowt.org>

5
This file is part of the GNU Fortran runtime library (libgfortran).
6 7 8 9

Libgfortran is free software; you can redistribute it and/or
modify it under the terms of the GNU General Public
License as published by the Free Software Foundation; either
10
version 3 of the License, or (at your option) any later version.
11 12 13 14 15 16

Libgfortran is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
GNU General Public License for more details.

17 18 19 20 21 22 23 24
Under Section 7 of GPL version 3, you are granted additional
permissions described in the GCC Runtime Library Exception, version
3.1, as published by the Free Software Foundation.

You should have received a copy of the GNU General Public License and
a copy of the GCC Runtime Library Exception along with this program;
see the files COPYING3 and COPYING.RUNTIME respectively.  If not, see
<http://www.gnu.org/licenses/>.  */
25

26
#include "libgfortran.h"
27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53
#include <stdlib.h>
#include <assert.h>


#if defined (HAVE_GFC_INTEGER_2) && defined (HAVE_GFC_INTEGER_2)


extern void minval_i2 (gfc_array_i2 * const restrict, 
	gfc_array_i2 * const restrict, const index_type * const restrict);
export_proto(minval_i2);

void
minval_i2 (gfc_array_i2 * const restrict retarray, 
	gfc_array_i2 * const restrict array, 
	const index_type * const restrict pdim)
{
  index_type count[GFC_MAX_DIMENSIONS];
  index_type extent[GFC_MAX_DIMENSIONS];
  index_type sstride[GFC_MAX_DIMENSIONS];
  index_type dstride[GFC_MAX_DIMENSIONS];
  const GFC_INTEGER_2 * restrict base;
  GFC_INTEGER_2 * restrict dest;
  index_type rank;
  index_type n;
  index_type len;
  index_type delta;
  index_type dim;
54
  int continue_loop;
55 56 57 58 59

  /* Make dim zero based to avoid confusion.  */
  dim = (*pdim) - 1;
  rank = GFC_DESCRIPTOR_RANK (array) - 1;

60
  len = GFC_DESCRIPTOR_EXTENT(array,dim);
61 62
  if (len < 0)
    len = 0;
63
  delta = GFC_DESCRIPTOR_STRIDE(array,dim);
64 65 66

  for (n = 0; n < dim; n++)
    {
67 68
      sstride[n] = GFC_DESCRIPTOR_STRIDE(array,n);
      extent[n] = GFC_DESCRIPTOR_EXTENT(array,n);
69 70 71 72 73 74

      if (extent[n] < 0)
	extent[n] = 0;
    }
  for (n = dim; n < rank; n++)
    {
75 76
      sstride[n] = GFC_DESCRIPTOR_STRIDE(array, n + 1);
      extent[n] = GFC_DESCRIPTOR_EXTENT(array, n + 1);
77 78 79 80 81

      if (extent[n] < 0)
	extent[n] = 0;
    }

82
  if (retarray->base_addr == NULL)
83
    {
84
      size_t alloc_size, str;
85 86

      for (n = 0; n < rank; n++)
87 88
	{
	  if (n == 0)
89
	    str = 1;
90 91
	  else
	    str = GFC_DESCRIPTOR_STRIDE(retarray,n-1) * extent[n-1];
92 93 94

	  GFC_DIMENSION_SET(retarray->dim[n], 0, extent[n] - 1, str);

95
	}
96 97 98 99

      retarray->offset = 0;
      retarray->dtype = (array->dtype & ~GFC_DTYPE_RANK_MASK) | rank;

100
      alloc_size = sizeof (GFC_INTEGER_2) * GFC_DESCRIPTOR_STRIDE(retarray,rank-1)
101 102
    		   * extent[rank-1];

103
      retarray->base_addr = xmalloc (alloc_size);
104 105 106
      if (alloc_size == 0)
	{
	  /* Make sure we have a zero-sized array.  */
107
	  GFC_DIMENSION_SET(retarray->dim[0], 0, -1, 1);
108
	  return;
109

110 111 112 113 114
	}
    }
  else
    {
      if (rank != GFC_DESCRIPTOR_RANK (retarray))
115
	runtime_error ("rank of return array incorrect in"
116 117 118
		       " MINVAL intrinsic: is %ld, should be %ld",
		       (long int) (GFC_DESCRIPTOR_RANK (retarray)),
		       (long int) rank);
119

120
      if (unlikely (compile_options.bounds_check))
Thomas Koenig committed
121 122
	bounds_ifunction_return ((array_t *) retarray, extent,
				 "return value", "MINVAL");
123 124 125 126 127
    }

  for (n = 0; n < rank; n++)
    {
      count[n] = 0;
128
      dstride[n] = GFC_DESCRIPTOR_STRIDE(retarray,n);
129
      if (extent[n] <= 0)
130
	return;
131 132
    }

133 134
  base = array->base_addr;
  dest = retarray->base_addr;
135

136 137
  continue_loop = 1;
  while (continue_loop)
138 139 140 141 142 143
    {
      const GFC_INTEGER_2 * restrict src;
      GFC_INTEGER_2 result;
      src = base;
      {

144 145 146 147 148 149
#if defined (GFC_INTEGER_2_INFINITY)
	result = GFC_INTEGER_2_INFINITY;
#else
	result = GFC_INTEGER_2_HUGE;
#endif
	if (len <= 0)
150 151 152 153 154 155
	  *dest = GFC_INTEGER_2_HUGE;
	else
	  {
	    for (n = 0; n < len; n++, src += delta)
	      {

156 157 158 159 160 161 162 163 164 165 166 167
#if defined (GFC_INTEGER_2_QUIET_NAN)
		if (*src <= result)
		  break;
	      }
	    if (unlikely (n >= len))
	      result = GFC_INTEGER_2_QUIET_NAN;
	    else for (; n < len; n++, src += delta)
	      {
#endif
		if (*src < result)
		  result = *src;
	      }
168
	    
169 170 171 172 173 174 175 176 177
	    *dest = result;
	  }
      }
      /* Advance to the next element.  */
      count[0]++;
      base += sstride[0];
      dest += dstride[0];
      n = 0;
      while (count[n] == extent[n])
178 179 180 181 182 183 184 185 186 187 188 189
	{
	  /* When we get to the end of a dimension, reset it and increment
	     the next dimension.  */
	  count[n] = 0;
	  /* We could precalculate these products, but this is a less
	     frequently used path so probably not worth it.  */
	  base -= sstride[n] * extent[n];
	  dest -= dstride[n] * extent[n];
	  n++;
	  if (n == rank)
	    {
	      /* Break out of the look.  */
190 191
	      continue_loop = 0;
	      break;
192 193 194 195 196 197 198 199
	    }
	  else
	    {
	      count[n]++;
	      base += sstride[n];
	      dest += dstride[n];
	    }
	}
200 201 202 203 204 205
    }
}


extern void mminval_i2 (gfc_array_i2 * const restrict, 
	gfc_array_i2 * const restrict, const index_type * const restrict,
206
	gfc_array_l1 * const restrict);
207 208 209 210 211 212
export_proto(mminval_i2);

void
mminval_i2 (gfc_array_i2 * const restrict retarray, 
	gfc_array_i2 * const restrict array, 
	const index_type * const restrict pdim, 
213
	gfc_array_l1 * const restrict mask)
214 215 216 217 218 219 220 221
{
  index_type count[GFC_MAX_DIMENSIONS];
  index_type extent[GFC_MAX_DIMENSIONS];
  index_type sstride[GFC_MAX_DIMENSIONS];
  index_type dstride[GFC_MAX_DIMENSIONS];
  index_type mstride[GFC_MAX_DIMENSIONS];
  GFC_INTEGER_2 * restrict dest;
  const GFC_INTEGER_2 * restrict base;
222
  const GFC_LOGICAL_1 * restrict mbase;
223 224 225 226 227 228
  int rank;
  int dim;
  index_type n;
  index_type len;
  index_type delta;
  index_type mdelta;
229
  int mask_kind;
230 231 232 233

  dim = (*pdim) - 1;
  rank = GFC_DESCRIPTOR_RANK (array) - 1;

234
  len = GFC_DESCRIPTOR_EXTENT(array,dim);
235 236
  if (len <= 0)
    return;
237

238
  mbase = mask->base_addr;
239 240 241 242 243 244 245 246 247 248 249 250

  mask_kind = GFC_DESCRIPTOR_SIZE (mask);

  if (mask_kind == 1 || mask_kind == 2 || mask_kind == 4 || mask_kind == 8
#ifdef HAVE_GFC_LOGICAL_16
      || mask_kind == 16
#endif
      )
    mbase = GFOR_POINTER_TO_L1 (mbase, mask_kind);
  else
    runtime_error ("Funny sized logical array");

251 252
  delta = GFC_DESCRIPTOR_STRIDE(array,dim);
  mdelta = GFC_DESCRIPTOR_STRIDE_BYTES(mask,dim);
253 254 255

  for (n = 0; n < dim; n++)
    {
256 257 258
      sstride[n] = GFC_DESCRIPTOR_STRIDE(array,n);
      mstride[n] = GFC_DESCRIPTOR_STRIDE_BYTES(mask,n);
      extent[n] = GFC_DESCRIPTOR_EXTENT(array,n);
259 260 261 262 263 264 265

      if (extent[n] < 0)
	extent[n] = 0;

    }
  for (n = dim; n < rank; n++)
    {
266 267 268
      sstride[n] = GFC_DESCRIPTOR_STRIDE(array,n + 1);
      mstride[n] = GFC_DESCRIPTOR_STRIDE_BYTES(mask, n + 1);
      extent[n] = GFC_DESCRIPTOR_EXTENT(array, n + 1);
269 270 271 272 273

      if (extent[n] < 0)
	extent[n] = 0;
    }

274
  if (retarray->base_addr == NULL)
275
    {
276
      size_t alloc_size, str;
277 278

      for (n = 0; n < rank; n++)
279 280 281 282 283
	{
	  if (n == 0)
	    str = 1;
	  else
	    str= GFC_DESCRIPTOR_STRIDE(retarray,n-1) * extent[n-1];
284 285 286

	  GFC_DIMENSION_SET(retarray->dim[n], 0, extent[n] - 1, str);

287
	}
288

289
      alloc_size = sizeof (GFC_INTEGER_2) * GFC_DESCRIPTOR_STRIDE(retarray,rank-1)
290 291 292 293 294 295 296 297
    		   * extent[rank-1];

      retarray->offset = 0;
      retarray->dtype = (array->dtype & ~GFC_DTYPE_RANK_MASK) | rank;

      if (alloc_size == 0)
	{
	  /* Make sure we have a zero-sized array.  */
298
	  GFC_DIMENSION_SET(retarray->dim[0], 0, -1, 1);
299 300 301
	  return;
	}
      else
302
	retarray->base_addr = xmalloc (alloc_size);
303 304 305 306 307

    }
  else
    {
      if (rank != GFC_DESCRIPTOR_RANK (retarray))
308 309
	runtime_error ("rank of return array incorrect in MINVAL intrinsic");

310
      if (unlikely (compile_options.bounds_check))
311
	{
Thomas Koenig committed
312 313 314 315
	  bounds_ifunction_return ((array_t *) retarray, extent,
				   "return value", "MINVAL");
	  bounds_equal_extents ((array_t *) mask, (array_t *) array,
	  			"MASK argument", "MINVAL");
316
	}
317 318 319 320 321
    }

  for (n = 0; n < rank; n++)
    {
      count[n] = 0;
322
      dstride[n] = GFC_DESCRIPTOR_STRIDE(retarray,n);
323
      if (extent[n] <= 0)
324
	return;
325 326
    }

327 328
  dest = retarray->base_addr;
  base = array->base_addr;
329 330 331 332

  while (base)
    {
      const GFC_INTEGER_2 * restrict src;
333
      const GFC_LOGICAL_1 * restrict msrc;
334 335 336 337 338
      GFC_INTEGER_2 result;
      src = base;
      msrc = mbase;
      {

339 340 341 342 343 344 345 346
#if defined (GFC_INTEGER_2_INFINITY)
	result = GFC_INTEGER_2_INFINITY;
#else
	result = GFC_INTEGER_2_HUGE;
#endif
#if defined (GFC_INTEGER_2_QUIET_NAN)
	int non_empty_p = 0;
#endif
347
	for (n = 0; n < len; n++, src += delta, msrc += mdelta)
348 349
	  {

350 351 352 353 354 355 356 357 358 359 360 361 362 363 364 365 366 367 368 369 370 371 372
#if defined (GFC_INTEGER_2_INFINITY) || defined (GFC_INTEGER_2_QUIET_NAN)
		if (*msrc)
		  {
#if defined (GFC_INTEGER_2_QUIET_NAN)
		    non_empty_p = 1;
		    if (*src <= result)
#endif
		      break;
		  }
	      }
	    if (unlikely (n >= len))
	      {
#if defined (GFC_INTEGER_2_QUIET_NAN)
		result = non_empty_p ? GFC_INTEGER_2_QUIET_NAN : GFC_INTEGER_2_HUGE;
#else
		result = GFC_INTEGER_2_HUGE;
#endif
	      }
	    else for (; n < len; n++, src += delta, msrc += mdelta)
	      {
#endif
		if (*msrc && *src < result)
		  result = *src;
373
	  }
374
	*dest = result;
375 376 377 378 379 380 381 382
      }
      /* Advance to the next element.  */
      count[0]++;
      base += sstride[0];
      mbase += mstride[0];
      dest += dstride[0];
      n = 0;
      while (count[n] == extent[n])
383 384 385 386 387 388 389 390 391 392 393 394 395 396 397 398 399 400 401 402 403 404 405 406
	{
	  /* When we get to the end of a dimension, reset it and increment
	     the next dimension.  */
	  count[n] = 0;
	  /* We could precalculate these products, but this is a less
	     frequently used path so probably not worth it.  */
	  base -= sstride[n] * extent[n];
	  mbase -= mstride[n] * extent[n];
	  dest -= dstride[n] * extent[n];
	  n++;
	  if (n == rank)
	    {
	      /* Break out of the look.  */
	      base = NULL;
	      break;
	    }
	  else
	    {
	      count[n]++;
	      base += sstride[n];
	      mbase += mstride[n];
	      dest += dstride[n];
	    }
	}
407 408 409 410 411 412 413 414 415 416 417 418 419 420 421
    }
}


extern void sminval_i2 (gfc_array_i2 * const restrict, 
	gfc_array_i2 * const restrict, const index_type * const restrict,
	GFC_LOGICAL_4 *);
export_proto(sminval_i2);

void
sminval_i2 (gfc_array_i2 * const restrict retarray, 
	gfc_array_i2 * const restrict array, 
	const index_type * const restrict pdim, 
	GFC_LOGICAL_4 * mask)
{
422 423 424 425
  index_type count[GFC_MAX_DIMENSIONS];
  index_type extent[GFC_MAX_DIMENSIONS];
  index_type dstride[GFC_MAX_DIMENSIONS];
  GFC_INTEGER_2 * restrict dest;
426 427
  index_type rank;
  index_type n;
428 429
  index_type dim;

430 431 432 433 434 435

  if (*mask)
    {
      minval_i2 (retarray, array, pdim);
      return;
    }
436 437 438 439 440 441
  /* Make dim zero based to avoid confusion.  */
  dim = (*pdim) - 1;
  rank = GFC_DESCRIPTOR_RANK (array) - 1;

  for (n = 0; n < dim; n++)
    {
442
      extent[n] = GFC_DESCRIPTOR_EXTENT(array,n);
443 444 445 446 447 448 449 450

      if (extent[n] <= 0)
	extent[n] = 0;
    }

  for (n = dim; n < rank; n++)
    {
      extent[n] =
451
	GFC_DESCRIPTOR_EXTENT(array,n + 1);
452 453

      if (extent[n] <= 0)
454
	extent[n] = 0;
455
    }
456

457
  if (retarray->base_addr == NULL)
458
    {
459
      size_t alloc_size, str;
460 461

      for (n = 0; n < rank; n++)
462 463 464 465 466
	{
	  if (n == 0)
	    str = 1;
	  else
	    str = GFC_DESCRIPTOR_STRIDE(retarray,n-1) * extent[n-1];
467 468 469

	  GFC_DIMENSION_SET(retarray->dim[n], 0, extent[n] - 1, str);

470
	}
471

472
      retarray->offset = 0;
473 474
      retarray->dtype = (array->dtype & ~GFC_DTYPE_RANK_MASK) | rank;

475
      alloc_size = sizeof (GFC_INTEGER_2) * GFC_DESCRIPTOR_STRIDE(retarray,rank-1)
476 477 478 479 480
    		   * extent[rank-1];

      if (alloc_size == 0)
	{
	  /* Make sure we have a zero-sized array.  */
481
	  GFC_DIMENSION_SET(retarray->dim[0], 0, -1, 1);
482 483 484
	  return;
	}
      else
485
	retarray->base_addr = xmalloc (alloc_size);
486 487 488
    }
  else
    {
489 490 491 492 493 494
      if (rank != GFC_DESCRIPTOR_RANK (retarray))
	runtime_error ("rank of return array incorrect in"
		       " MINVAL intrinsic: is %ld, should be %ld",
		       (long int) (GFC_DESCRIPTOR_RANK (retarray)),
		       (long int) rank);

495
      if (unlikely (compile_options.bounds_check))
496
	{
497 498 499
	  for (n=0; n < rank; n++)
	    {
	      index_type ret_extent;
500

501
	      ret_extent = GFC_DESCRIPTOR_EXTENT(retarray,n);
502 503 504 505 506 507
	      if (extent[n] != ret_extent)
		runtime_error ("Incorrect extent in return value of"
			       " MINVAL intrinsic in dimension %ld:"
			       " is %ld, should be %ld", (long int) n + 1,
			       (long int) ret_extent, (long int) extent[n]);
	    }
508 509
	}
    }
510

511 512 513
  for (n = 0; n < rank; n++)
    {
      count[n] = 0;
514
      dstride[n] = GFC_DESCRIPTOR_STRIDE(retarray,n);
515 516
    }

517
  dest = retarray->base_addr;
518 519 520 521 522 523 524 525

  while(1)
    {
      *dest = GFC_INTEGER_2_HUGE;
      count[0]++;
      dest += dstride[0];
      n = 0;
      while (count[n] == extent[n])
526
	{
527
	  /* When we get to the end of a dimension, reset it and increment
528 529 530 531 532 533 534
	     the next dimension.  */
	  count[n] = 0;
	  /* We could precalculate these products, but this is a less
	     frequently used path so probably not worth it.  */
	  dest -= dstride[n] * extent[n];
	  n++;
	  if (n == rank)
535
	    return;
536 537 538 539 540
	  else
	    {
	      count[n]++;
	      dest += dstride[n];
	    }
541 542
      	}
    }
543 544 545
}

#endif