1 /* Implementation of the MAXVAL intrinsic
2 Copyright (C) 2002-2013 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 <stdlib.h>
28 #include <assert.h>
29
30
31 #if defined (HAVE_GFC_REAL_8) && defined (HAVE_GFC_REAL_8)
32
33
34 extern void maxval_r8 (gfc_array_r8 * const restrict,
35 gfc_array_r8 * const restrict, const index_type * const restrict);
36 export_proto(maxval_r8);
37
38 void
maxval_r8(gfc_array_r8 * const restrict retarray,gfc_array_r8 * const restrict array,const index_type * const restrict pdim)39 maxval_r8 (gfc_array_r8 * const restrict retarray,
40 gfc_array_r8 * const restrict array,
41 const index_type * const restrict pdim)
42 {
43 index_type count[GFC_MAX_DIMENSIONS];
44 index_type extent[GFC_MAX_DIMENSIONS];
45 index_type sstride[GFC_MAX_DIMENSIONS];
46 index_type dstride[GFC_MAX_DIMENSIONS];
47 const GFC_REAL_8 * restrict base;
48 GFC_REAL_8 * restrict dest;
49 index_type rank;
50 index_type n;
51 index_type len;
52 index_type delta;
53 index_type dim;
54 int continue_loop;
55
56 /* Make dim zero based to avoid confusion. */
57 dim = (*pdim) - 1;
58 rank = GFC_DESCRIPTOR_RANK (array) - 1;
59
60 len = GFC_DESCRIPTOR_EXTENT(array,dim);
61 if (len < 0)
62 len = 0;
63 delta = GFC_DESCRIPTOR_STRIDE(array,dim);
64
65 for (n = 0; n < dim; n++)
66 {
67 sstride[n] = GFC_DESCRIPTOR_STRIDE(array,n);
68 extent[n] = GFC_DESCRIPTOR_EXTENT(array,n);
69
70 if (extent[n] < 0)
71 extent[n] = 0;
72 }
73 for (n = dim; n < rank; n++)
74 {
75 sstride[n] = GFC_DESCRIPTOR_STRIDE(array, n + 1);
76 extent[n] = GFC_DESCRIPTOR_EXTENT(array, n + 1);
77
78 if (extent[n] < 0)
79 extent[n] = 0;
80 }
81
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 else
91 str = GFC_DESCRIPTOR_STRIDE(retarray,n-1) * extent[n-1];
92
93 GFC_DIMENSION_SET(retarray->dim[n], 0, extent[n] - 1, str);
94
95 }
96
97 retarray->offset = 0;
98 retarray->dtype = (array->dtype & ~GFC_DTYPE_RANK_MASK) | rank;
99
100 alloc_size = GFC_DESCRIPTOR_STRIDE(retarray,rank-1) * extent[rank-1];
101
102 retarray->base_addr = xmallocarray (alloc_size, sizeof (GFC_REAL_8));
103 if (alloc_size == 0)
104 {
105 /* Make sure we have a zero-sized array. */
106 GFC_DIMENSION_SET(retarray->dim[0], 0, -1, 1);
107 return;
108
109 }
110 }
111 else
112 {
113 if (rank != GFC_DESCRIPTOR_RANK (retarray))
114 runtime_error ("rank of return array incorrect in"
115 " MAXVAL intrinsic: is %ld, should be %ld",
116 (long int) (GFC_DESCRIPTOR_RANK (retarray)),
117 (long int) rank);
118
119 if (unlikely (compile_options.bounds_check))
120 bounds_ifunction_return ((array_t *) retarray, extent,
121 "return value", "MAXVAL");
122 }
123
124 for (n = 0; n < rank; n++)
125 {
126 count[n] = 0;
127 dstride[n] = GFC_DESCRIPTOR_STRIDE(retarray,n);
128 if (extent[n] <= 0)
129 return;
130 }
131
132 base = array->base_addr;
133 dest = retarray->base_addr;
134
135 continue_loop = 1;
136 while (continue_loop)
137 {
138 const GFC_REAL_8 * restrict src;
139 GFC_REAL_8 result;
140 src = base;
141 {
142
143 #if defined (GFC_REAL_8_INFINITY)
144 result = -GFC_REAL_8_INFINITY;
145 #else
146 result = -GFC_REAL_8_HUGE;
147 #endif
148 if (len <= 0)
149 *dest = -GFC_REAL_8_HUGE;
150 else
151 {
152 for (n = 0; n < len; n++, src += delta)
153 {
154
155 #if defined (GFC_REAL_8_QUIET_NAN)
156 if (*src >= result)
157 break;
158 }
159 if (unlikely (n >= len))
160 result = GFC_REAL_8_QUIET_NAN;
161 else for (; n < len; n++, src += delta)
162 {
163 #endif
164 if (*src > result)
165 result = *src;
166 }
167
168 *dest = result;
169 }
170 }
171 /* Advance to the next element. */
172 count[0]++;
173 base += sstride[0];
174 dest += dstride[0];
175 n = 0;
176 while (count[n] == extent[n])
177 {
178 /* When we get to the end of a dimension, reset it and increment
179 the next dimension. */
180 count[n] = 0;
181 /* We could precalculate these products, but this is a less
182 frequently used path so probably not worth it. */
183 base -= sstride[n] * extent[n];
184 dest -= dstride[n] * extent[n];
185 n++;
186 if (n == rank)
187 {
188 /* Break out of the look. */
189 continue_loop = 0;
190 break;
191 }
192 else
193 {
194 count[n]++;
195 base += sstride[n];
196 dest += dstride[n];
197 }
198 }
199 }
200 }
201
202
203 extern void mmaxval_r8 (gfc_array_r8 * const restrict,
204 gfc_array_r8 * const restrict, const index_type * const restrict,
205 gfc_array_l1 * const restrict);
206 export_proto(mmaxval_r8);
207
208 void
mmaxval_r8(gfc_array_r8 * const restrict retarray,gfc_array_r8 * const restrict array,const index_type * const restrict pdim,gfc_array_l1 * const restrict mask)209 mmaxval_r8 (gfc_array_r8 * const restrict retarray,
210 gfc_array_r8 * const restrict array,
211 const index_type * const restrict pdim,
212 gfc_array_l1 * const restrict mask)
213 {
214 index_type count[GFC_MAX_DIMENSIONS];
215 index_type extent[GFC_MAX_DIMENSIONS];
216 index_type sstride[GFC_MAX_DIMENSIONS];
217 index_type dstride[GFC_MAX_DIMENSIONS];
218 index_type mstride[GFC_MAX_DIMENSIONS];
219 GFC_REAL_8 * restrict dest;
220 const GFC_REAL_8 * restrict base;
221 const GFC_LOGICAL_1 * restrict mbase;
222 int rank;
223 int dim;
224 index_type n;
225 index_type len;
226 index_type delta;
227 index_type mdelta;
228 int mask_kind;
229
230 dim = (*pdim) - 1;
231 rank = GFC_DESCRIPTOR_RANK (array) - 1;
232
233 len = GFC_DESCRIPTOR_EXTENT(array,dim);
234 if (len <= 0)
235 return;
236
237 mbase = mask->base_addr;
238
239 mask_kind = GFC_DESCRIPTOR_SIZE (mask);
240
241 if (mask_kind == 1 || mask_kind == 2 || mask_kind == 4 || mask_kind == 8
242 #ifdef HAVE_GFC_LOGICAL_16
243 || mask_kind == 16
244 #endif
245 )
246 mbase = GFOR_POINTER_TO_L1 (mbase, mask_kind);
247 else
248 runtime_error ("Funny sized logical array");
249
250 delta = GFC_DESCRIPTOR_STRIDE(array,dim);
251 mdelta = GFC_DESCRIPTOR_STRIDE_BYTES(mask,dim);
252
253 for (n = 0; n < dim; n++)
254 {
255 sstride[n] = GFC_DESCRIPTOR_STRIDE(array,n);
256 mstride[n] = GFC_DESCRIPTOR_STRIDE_BYTES(mask,n);
257 extent[n] = GFC_DESCRIPTOR_EXTENT(array,n);
258
259 if (extent[n] < 0)
260 extent[n] = 0;
261
262 }
263 for (n = dim; n < rank; n++)
264 {
265 sstride[n] = GFC_DESCRIPTOR_STRIDE(array,n + 1);
266 mstride[n] = GFC_DESCRIPTOR_STRIDE_BYTES(mask, n + 1);
267 extent[n] = GFC_DESCRIPTOR_EXTENT(array, n + 1);
268
269 if (extent[n] < 0)
270 extent[n] = 0;
271 }
272
273 if (retarray->base_addr == NULL)
274 {
275 size_t alloc_size, str;
276
277 for (n = 0; n < rank; n++)
278 {
279 if (n == 0)
280 str = 1;
281 else
282 str= GFC_DESCRIPTOR_STRIDE(retarray,n-1) * extent[n-1];
283
284 GFC_DIMENSION_SET(retarray->dim[n], 0, extent[n] - 1, str);
285
286 }
287
288 alloc_size = GFC_DESCRIPTOR_STRIDE(retarray,rank-1) * extent[rank-1];
289
290 retarray->offset = 0;
291 retarray->dtype = (array->dtype & ~GFC_DTYPE_RANK_MASK) | rank;
292
293 if (alloc_size == 0)
294 {
295 /* Make sure we have a zero-sized array. */
296 GFC_DIMENSION_SET(retarray->dim[0], 0, -1, 1);
297 return;
298 }
299 else
300 retarray->base_addr = xmallocarray (alloc_size, sizeof (GFC_REAL_8));
301
302 }
303 else
304 {
305 if (rank != GFC_DESCRIPTOR_RANK (retarray))
306 runtime_error ("rank of return array incorrect in MAXVAL intrinsic");
307
308 if (unlikely (compile_options.bounds_check))
309 {
310 bounds_ifunction_return ((array_t *) retarray, extent,
311 "return value", "MAXVAL");
312 bounds_equal_extents ((array_t *) mask, (array_t *) array,
313 "MASK argument", "MAXVAL");
314 }
315 }
316
317 for (n = 0; n < rank; n++)
318 {
319 count[n] = 0;
320 dstride[n] = GFC_DESCRIPTOR_STRIDE(retarray,n);
321 if (extent[n] <= 0)
322 return;
323 }
324
325 dest = retarray->base_addr;
326 base = array->base_addr;
327
328 while (base)
329 {
330 const GFC_REAL_8 * restrict src;
331 const GFC_LOGICAL_1 * restrict msrc;
332 GFC_REAL_8 result;
333 src = base;
334 msrc = mbase;
335 {
336
337 #if defined (GFC_REAL_8_INFINITY)
338 result = -GFC_REAL_8_INFINITY;
339 #else
340 result = -GFC_REAL_8_HUGE;
341 #endif
342 #if defined (GFC_REAL_8_QUIET_NAN)
343 int non_empty_p = 0;
344 #endif
345 for (n = 0; n < len; n++, src += delta, msrc += mdelta)
346 {
347
348 #if defined (GFC_REAL_8_INFINITY) || defined (GFC_REAL_8_QUIET_NAN)
349 if (*msrc)
350 {
351 #if defined (GFC_REAL_8_QUIET_NAN)
352 non_empty_p = 1;
353 if (*src >= result)
354 #endif
355 break;
356 }
357 }
358 if (unlikely (n >= len))
359 {
360 #if defined (GFC_REAL_8_QUIET_NAN)
361 result = non_empty_p ? GFC_REAL_8_QUIET_NAN : -GFC_REAL_8_HUGE;
362 #else
363 result = -GFC_REAL_8_HUGE;
364 #endif
365 }
366 else for (; n < len; n++, src += delta, msrc += mdelta)
367 {
368 #endif
369 if (*msrc && *src > result)
370 result = *src;
371 }
372 *dest = result;
373 }
374 /* Advance to the next element. */
375 count[0]++;
376 base += sstride[0];
377 mbase += mstride[0];
378 dest += dstride[0];
379 n = 0;
380 while (count[n] == extent[n])
381 {
382 /* When we get to the end of a dimension, reset it and increment
383 the next dimension. */
384 count[n] = 0;
385 /* We could precalculate these products, but this is a less
386 frequently used path so probably not worth it. */
387 base -= sstride[n] * extent[n];
388 mbase -= mstride[n] * extent[n];
389 dest -= dstride[n] * extent[n];
390 n++;
391 if (n == rank)
392 {
393 /* Break out of the look. */
394 base = NULL;
395 break;
396 }
397 else
398 {
399 count[n]++;
400 base += sstride[n];
401 mbase += mstride[n];
402 dest += dstride[n];
403 }
404 }
405 }
406 }
407
408
409 extern void smaxval_r8 (gfc_array_r8 * const restrict,
410 gfc_array_r8 * const restrict, const index_type * const restrict,
411 GFC_LOGICAL_4 *);
412 export_proto(smaxval_r8);
413
414 void
smaxval_r8(gfc_array_r8 * const restrict retarray,gfc_array_r8 * const restrict array,const index_type * const restrict pdim,GFC_LOGICAL_4 * mask)415 smaxval_r8 (gfc_array_r8 * const restrict retarray,
416 gfc_array_r8 * const restrict array,
417 const index_type * const restrict pdim,
418 GFC_LOGICAL_4 * mask)
419 {
420 index_type count[GFC_MAX_DIMENSIONS];
421 index_type extent[GFC_MAX_DIMENSIONS];
422 index_type dstride[GFC_MAX_DIMENSIONS];
423 GFC_REAL_8 * restrict dest;
424 index_type rank;
425 index_type n;
426 index_type dim;
427
428
429 if (*mask)
430 {
431 maxval_r8 (retarray, array, pdim);
432 return;
433 }
434 /* Make dim zero based to avoid confusion. */
435 dim = (*pdim) - 1;
436 rank = GFC_DESCRIPTOR_RANK (array) - 1;
437
438 for (n = 0; n < dim; n++)
439 {
440 extent[n] = GFC_DESCRIPTOR_EXTENT(array,n);
441
442 if (extent[n] <= 0)
443 extent[n] = 0;
444 }
445
446 for (n = dim; n < rank; n++)
447 {
448 extent[n] =
449 GFC_DESCRIPTOR_EXTENT(array,n + 1);
450
451 if (extent[n] <= 0)
452 extent[n] = 0;
453 }
454
455 if (retarray->base_addr == NULL)
456 {
457 size_t alloc_size, str;
458
459 for (n = 0; n < rank; n++)
460 {
461 if (n == 0)
462 str = 1;
463 else
464 str = GFC_DESCRIPTOR_STRIDE(retarray,n-1) * extent[n-1];
465
466 GFC_DIMENSION_SET(retarray->dim[n], 0, extent[n] - 1, str);
467
468 }
469
470 retarray->offset = 0;
471 retarray->dtype = (array->dtype & ~GFC_DTYPE_RANK_MASK) | rank;
472
473 alloc_size = GFC_DESCRIPTOR_STRIDE(retarray,rank-1) * extent[rank-1];
474
475 if (alloc_size == 0)
476 {
477 /* Make sure we have a zero-sized array. */
478 GFC_DIMENSION_SET(retarray->dim[0], 0, -1, 1);
479 return;
480 }
481 else
482 retarray->base_addr = xmallocarray (alloc_size, sizeof (GFC_REAL_8));
483 }
484 else
485 {
486 if (rank != GFC_DESCRIPTOR_RANK (retarray))
487 runtime_error ("rank of return array incorrect in"
488 " MAXVAL intrinsic: is %ld, should be %ld",
489 (long int) (GFC_DESCRIPTOR_RANK (retarray)),
490 (long int) rank);
491
492 if (unlikely (compile_options.bounds_check))
493 {
494 for (n=0; n < rank; n++)
495 {
496 index_type ret_extent;
497
498 ret_extent = GFC_DESCRIPTOR_EXTENT(retarray,n);
499 if (extent[n] != ret_extent)
500 runtime_error ("Incorrect extent in return value of"
501 " MAXVAL intrinsic in dimension %ld:"
502 " is %ld, should be %ld", (long int) n + 1,
503 (long int) ret_extent, (long int) extent[n]);
504 }
505 }
506 }
507
508 for (n = 0; n < rank; n++)
509 {
510 count[n] = 0;
511 dstride[n] = GFC_DESCRIPTOR_STRIDE(retarray,n);
512 }
513
514 dest = retarray->base_addr;
515
516 while(1)
517 {
518 *dest = -GFC_REAL_8_HUGE;
519 count[0]++;
520 dest += dstride[0];
521 n = 0;
522 while (count[n] == extent[n])
523 {
524 /* When we get to the end of a dimension, reset it and increment
525 the next dimension. */
526 count[n] = 0;
527 /* We could precalculate these products, but this is a less
528 frequently used path so probably not worth it. */
529 dest -= dstride[n] * extent[n];
530 n++;
531 if (n == rank)
532 return;
533 else
534 {
535 count[n]++;
536 dest += dstride[n];
537 }
538 }
539 }
540 }
541
542 #endif
543