1 /* Implementation of the PRODUCT intrinsic
2 Copyright 2002 Free Software Foundation, Inc.
3 Contributed by Paul Brook <paul@nowt.org>
5 This file is part of the GNU Fortran 95 runtime library (libgfortran).
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 2 of the License, or (at your option) any later version.
12 In addition to the permissions in the GNU General Public License, the
13 Free Software Foundation gives you unlimited permission to link the
14 compiled version of this file into combinations with other programs,
15 and to distribute those combinations without any restriction coming
16 from the use of this file. (The General Public License restrictions
17 do apply in other respects; for example, they cover modification of
18 the file, and distribution when not linked into a combine
21 Libgfortran is distributed in the hope that it will be useful,
22 but WITHOUT ANY WARRANTY; without even the implied warranty of
23 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
24 GNU General Public License for more details.
26 You should have received a copy of the GNU General Public
27 License along with libgfortran; see the file COPYING. If not,
28 write to the Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor,
29 Boston, MA 02110-1301, USA. */
34 #include "libgfortran.h"
37 #if defined (HAVE_GFC_REAL_16) && defined (HAVE_GFC_REAL_16)
40 extern void product_r16 (gfc_array_r16
* const restrict
,
41 gfc_array_r16
* const restrict
, const index_type
* const restrict
);
42 export_proto(product_r16
);
45 product_r16 (gfc_array_r16
* const restrict retarray
,
46 gfc_array_r16
* const restrict array
,
47 const index_type
* const restrict pdim
)
49 index_type count
[GFC_MAX_DIMENSIONS
];
50 index_type extent
[GFC_MAX_DIMENSIONS
];
51 index_type sstride
[GFC_MAX_DIMENSIONS
];
52 index_type dstride
[GFC_MAX_DIMENSIONS
];
53 const GFC_REAL_16
* restrict base
;
54 GFC_REAL_16
* restrict dest
;
61 /* Make dim zero based to avoid confusion. */
63 rank
= GFC_DESCRIPTOR_RANK (array
) - 1;
65 len
= array
->dim
[dim
].ubound
+ 1 - array
->dim
[dim
].lbound
;
66 delta
= array
->dim
[dim
].stride
;
68 for (n
= 0; n
< dim
; n
++)
70 sstride
[n
] = array
->dim
[n
].stride
;
71 extent
[n
] = array
->dim
[n
].ubound
+ 1 - array
->dim
[n
].lbound
;
76 for (n
= dim
; n
< rank
; n
++)
78 sstride
[n
] = array
->dim
[n
+ 1].stride
;
80 array
->dim
[n
+ 1].ubound
+ 1 - array
->dim
[n
+ 1].lbound
;
86 if (retarray
->data
== NULL
)
90 for (n
= 0; n
< rank
; n
++)
92 retarray
->dim
[n
].lbound
= 0;
93 retarray
->dim
[n
].ubound
= extent
[n
]-1;
95 retarray
->dim
[n
].stride
= 1;
97 retarray
->dim
[n
].stride
= retarray
->dim
[n
-1].stride
* extent
[n
-1];
100 retarray
->offset
= 0;
101 retarray
->dtype
= (array
->dtype
& ~GFC_DTYPE_RANK_MASK
) | rank
;
103 alloc_size
= sizeof (GFC_REAL_16
) * retarray
->dim
[rank
-1].stride
108 /* Make sure we have a zero-sized array. */
109 retarray
->dim
[0].lbound
= 0;
110 retarray
->dim
[0].ubound
= -1;
114 retarray
->data
= internal_malloc_size (alloc_size
);
118 if (rank
!= GFC_DESCRIPTOR_RANK (retarray
))
119 runtime_error ("rank of return array incorrect");
122 for (n
= 0; n
< rank
; n
++)
125 dstride
[n
] = retarray
->dim
[n
].stride
;
131 dest
= retarray
->data
;
135 const GFC_REAL_16
* restrict src
;
145 for (n
= 0; n
< len
; n
++, src
+= delta
)
153 /* Advance to the next element. */
158 while (count
[n
] == extent
[n
])
160 /* When we get to the end of a dimension, reset it and increment
161 the next dimension. */
163 /* We could precalculate these products, but this is a less
164 frequently used path so probably not worth it. */
165 base
-= sstride
[n
] * extent
[n
];
166 dest
-= dstride
[n
] * extent
[n
];
170 /* Break out of the look. */
185 extern void mproduct_r16 (gfc_array_r16
* const restrict
,
186 gfc_array_r16
* const restrict
, const index_type
* const restrict
,
187 gfc_array_l1
* const restrict
);
188 export_proto(mproduct_r16
);
191 mproduct_r16 (gfc_array_r16
* const restrict retarray
,
192 gfc_array_r16
* const restrict array
,
193 const index_type
* const restrict pdim
,
194 gfc_array_l1
* const restrict mask
)
196 index_type count
[GFC_MAX_DIMENSIONS
];
197 index_type extent
[GFC_MAX_DIMENSIONS
];
198 index_type sstride
[GFC_MAX_DIMENSIONS
];
199 index_type dstride
[GFC_MAX_DIMENSIONS
];
200 index_type mstride
[GFC_MAX_DIMENSIONS
];
201 GFC_REAL_16
* restrict dest
;
202 const GFC_REAL_16
* restrict base
;
203 const GFC_LOGICAL_1
* restrict mbase
;
213 rank
= GFC_DESCRIPTOR_RANK (array
) - 1;
215 len
= array
->dim
[dim
].ubound
+ 1 - array
->dim
[dim
].lbound
;
221 mask_kind
= GFC_DESCRIPTOR_SIZE (mask
);
223 if (mask_kind
== 1 || mask_kind
== 2 || mask_kind
== 4 || mask_kind
== 8
224 #ifdef HAVE_GFC_LOGICAL_16
228 mbase
= GFOR_POINTER_TO_L1 (mbase
, mask_kind
);
230 runtime_error ("Funny sized logical array");
232 delta
= array
->dim
[dim
].stride
;
233 mdelta
= mask
->dim
[dim
].stride
* mask_kind
;
235 for (n
= 0; n
< dim
; n
++)
237 sstride
[n
] = array
->dim
[n
].stride
;
238 mstride
[n
] = mask
->dim
[n
].stride
* mask_kind
;
239 extent
[n
] = array
->dim
[n
].ubound
+ 1 - array
->dim
[n
].lbound
;
245 for (n
= dim
; n
< rank
; n
++)
247 sstride
[n
] = array
->dim
[n
+ 1].stride
;
248 mstride
[n
] = mask
->dim
[n
+ 1].stride
* mask_kind
;
250 array
->dim
[n
+ 1].ubound
+ 1 - array
->dim
[n
+ 1].lbound
;
256 if (retarray
->data
== NULL
)
260 for (n
= 0; n
< rank
; n
++)
262 retarray
->dim
[n
].lbound
= 0;
263 retarray
->dim
[n
].ubound
= extent
[n
]-1;
265 retarray
->dim
[n
].stride
= 1;
267 retarray
->dim
[n
].stride
= retarray
->dim
[n
-1].stride
* extent
[n
-1];
270 alloc_size
= sizeof (GFC_REAL_16
) * retarray
->dim
[rank
-1].stride
273 retarray
->offset
= 0;
274 retarray
->dtype
= (array
->dtype
& ~GFC_DTYPE_RANK_MASK
) | rank
;
278 /* Make sure we have a zero-sized array. */
279 retarray
->dim
[0].lbound
= 0;
280 retarray
->dim
[0].ubound
= -1;
284 retarray
->data
= internal_malloc_size (alloc_size
);
289 if (rank
!= GFC_DESCRIPTOR_RANK (retarray
))
290 runtime_error ("rank of return array incorrect");
293 for (n
= 0; n
< rank
; n
++)
296 dstride
[n
] = retarray
->dim
[n
].stride
;
301 dest
= retarray
->data
;
306 const GFC_REAL_16
* restrict src
;
307 const GFC_LOGICAL_1
* restrict msrc
;
318 for (n
= 0; n
< len
; n
++, src
+= delta
, msrc
+= mdelta
)
327 /* Advance to the next element. */
333 while (count
[n
] == extent
[n
])
335 /* When we get to the end of a dimension, reset it and increment
336 the next dimension. */
338 /* We could precalculate these products, but this is a less
339 frequently used path so probably not worth it. */
340 base
-= sstride
[n
] * extent
[n
];
341 mbase
-= mstride
[n
] * extent
[n
];
342 dest
-= dstride
[n
] * extent
[n
];
346 /* Break out of the look. */
362 extern void sproduct_r16 (gfc_array_r16
* const restrict
,
363 gfc_array_r16
* const restrict
, const index_type
* const restrict
,
365 export_proto(sproduct_r16
);
368 sproduct_r16 (gfc_array_r16
* const restrict retarray
,
369 gfc_array_r16
* const restrict array
,
370 const index_type
* const restrict pdim
,
371 GFC_LOGICAL_4
* mask
)
380 product_r16 (retarray
, array
, pdim
);
383 rank
= GFC_DESCRIPTOR_RANK (array
);
385 runtime_error ("Rank of array needs to be > 0");
387 if (retarray
->data
== NULL
)
389 retarray
->dim
[0].lbound
= 0;
390 retarray
->dim
[0].ubound
= rank
-1;
391 retarray
->dim
[0].stride
= 1;
392 retarray
->dtype
= (retarray
->dtype
& ~GFC_DTYPE_RANK_MASK
) | 1;
393 retarray
->offset
= 0;
394 retarray
->data
= internal_malloc_size (sizeof (GFC_REAL_16
) * rank
);
398 if (GFC_DESCRIPTOR_RANK (retarray
) != 1)
399 runtime_error ("rank of return array does not equal 1");
401 if (retarray
->dim
[0].ubound
+ 1 - retarray
->dim
[0].lbound
!= rank
)
402 runtime_error ("dimension of return array incorrect");
405 dstride
= retarray
->dim
[0].stride
;
406 dest
= retarray
->data
;
408 for (n
= 0; n
< rank
; n
++)
409 dest
[n
* dstride
] = 1 ;