]> git.ipfire.org Git - thirdparty/gcc.git/blob - libgfortran/generated/reshape_r4.c
re PR fortran/38135 (RESHAPE gives wrong result)
[thirdparty/gcc.git] / libgfortran / generated / reshape_r4.c
1 /* Implementation of the RESHAPE
2 Copyright 2002, 2006, 2007 Free Software Foundation, Inc.
3 Contributed by Paul Brook <paul@nowt.org>
4
5 This file is part of the GNU Fortran 95 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 2 of the License, or (at your option) any later version.
11
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
19 executable.)
20
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.
25
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. */
30
31 #include "libgfortran.h"
32 #include <stdlib.h>
33 #include <assert.h>
34
35
36 #if defined (HAVE_GFC_REAL_4)
37
38 typedef GFC_ARRAY_DESCRIPTOR(1, index_type) shape_type;
39
40
41 extern void reshape_r4 (gfc_array_r4 * const restrict,
42 gfc_array_r4 * const restrict,
43 shape_type * const restrict,
44 gfc_array_r4 * const restrict,
45 shape_type * const restrict);
46 export_proto(reshape_r4);
47
48 void
49 reshape_r4 (gfc_array_r4 * const restrict ret,
50 gfc_array_r4 * const restrict source,
51 shape_type * const restrict shape,
52 gfc_array_r4 * const restrict pad,
53 shape_type * const restrict order)
54 {
55 /* r.* indicates the return array. */
56 index_type rcount[GFC_MAX_DIMENSIONS];
57 index_type rextent[GFC_MAX_DIMENSIONS];
58 index_type rstride[GFC_MAX_DIMENSIONS];
59 index_type rstride0;
60 index_type rdim;
61 index_type rsize;
62 index_type rs;
63 index_type rex;
64 GFC_REAL_4 *rptr;
65 /* s.* indicates the source array. */
66 index_type scount[GFC_MAX_DIMENSIONS];
67 index_type sextent[GFC_MAX_DIMENSIONS];
68 index_type sstride[GFC_MAX_DIMENSIONS];
69 index_type sstride0;
70 index_type sdim;
71 index_type ssize;
72 const GFC_REAL_4 *sptr;
73 /* p.* indicates the pad array. */
74 index_type pcount[GFC_MAX_DIMENSIONS];
75 index_type pextent[GFC_MAX_DIMENSIONS];
76 index_type pstride[GFC_MAX_DIMENSIONS];
77 index_type pdim;
78 index_type psize;
79 const GFC_REAL_4 *pptr;
80
81 const GFC_REAL_4 *src;
82 int n;
83 int dim;
84 int sempty, pempty, shape_empty;
85 index_type shape_data[GFC_MAX_DIMENSIONS];
86
87 rdim = shape->dim[0].ubound - shape->dim[0].lbound + 1;
88 if (rdim != GFC_DESCRIPTOR_RANK(ret))
89 runtime_error("rank of return array incorrect in RESHAPE intrinsic");
90
91 shape_empty = 0;
92
93 for (n = 0; n < rdim; n++)
94 {
95 shape_data[n] = shape->data[n * shape->dim[0].stride];
96 if (shape_data[n] <= 0)
97 {
98 shape_data[n] = 0;
99 shape_empty = 1;
100 }
101 }
102
103 if (ret->data == NULL)
104 {
105 rs = 1;
106 for (n = 0; n < rdim; n++)
107 {
108 ret->dim[n].lbound = 0;
109 rex = shape_data[n];
110 ret->dim[n].ubound = rex - 1;
111 ret->dim[n].stride = rs;
112 rs *= rex;
113 }
114 ret->offset = 0;
115 ret->data = internal_malloc_size ( rs * sizeof (GFC_REAL_4));
116 ret->dtype = (source->dtype & ~GFC_DTYPE_RANK_MASK) | rdim;
117 }
118
119 if (shape_empty)
120 return;
121
122 if (unlikely (compile_options.bounds_check))
123 {
124 index_type ret_extent, source_extent;
125
126 rs = 1;
127 for (n = 0; n < rdim; n++)
128 {
129 rs *= shape_data[n];
130 ret_extent = ret->dim[n].ubound + 1 - ret->dim[n].lbound;
131 if (ret_extent != shape_data[n])
132 runtime_error("Incorrect extent in return value of RESHAPE"
133 " intrinsic in dimension %ld: is %ld,"
134 " should be %ld", (long int) n+1,
135 (long int) ret_extent, (long int) shape_data[n]);
136 }
137
138 source_extent = source->dim[0].ubound + 1 - source->dim[0].lbound;
139
140 if (rs < source_extent || (rs > source_extent && !pad))
141 runtime_error("Incorrect size in SOURCE argument to RESHAPE"
142 " intrinsic: is %ld, should be %ld",
143 (long int) source_extent, (long int) rs);
144
145 if (order)
146 {
147 int seen[GFC_MAX_DIMENSIONS];
148 index_type v;
149
150 for (n = 0; n < rdim; n++)
151 seen[n] = 0;
152
153 for (n = 0; n < rdim; n++)
154 {
155 v = order->data[n * order->dim[0].stride] - 1;
156
157 if (v < 0 || v >= rdim)
158 runtime_error("Value %ld out of range in ORDER argument"
159 " to RESHAPE intrinsic", (long int) v + 1);
160
161 if (seen[v] != 0)
162 runtime_error("Duplicate value %ld in ORDER argument to"
163 " RESHAPE intrinsic", (long int) v + 1);
164
165 seen[v] = 1;
166 }
167 }
168 }
169
170 rsize = 1;
171 for (n = 0; n < rdim; n++)
172 {
173 if (order)
174 dim = order->data[n * order->dim[0].stride] - 1;
175 else
176 dim = n;
177
178 rcount[n] = 0;
179 rstride[n] = ret->dim[dim].stride;
180 rextent[n] = ret->dim[dim].ubound + 1 - ret->dim[dim].lbound;
181 if (rextent[n] < 0)
182 rextent[n] = 0;
183
184 if (rextent[n] != shape_data[dim])
185 runtime_error ("shape and target do not conform");
186
187 if (rsize == rstride[n])
188 rsize *= rextent[n];
189 else
190 rsize = 0;
191 if (rextent[n] <= 0)
192 return;
193 }
194
195 sdim = GFC_DESCRIPTOR_RANK (source);
196 ssize = 1;
197 sempty = 0;
198 for (n = 0; n < sdim; n++)
199 {
200 scount[n] = 0;
201 sstride[n] = source->dim[n].stride;
202 sextent[n] = source->dim[n].ubound + 1 - source->dim[n].lbound;
203 if (sextent[n] <= 0)
204 {
205 sempty = 1;
206 sextent[n] = 0;
207 }
208
209 if (ssize == sstride[n])
210 ssize *= sextent[n];
211 else
212 ssize = 0;
213 }
214
215 if (pad)
216 {
217 pdim = GFC_DESCRIPTOR_RANK (pad);
218 psize = 1;
219 pempty = 0;
220 for (n = 0; n < pdim; n++)
221 {
222 pcount[n] = 0;
223 pstride[n] = pad->dim[n].stride;
224 pextent[n] = pad->dim[n].ubound + 1 - pad->dim[n].lbound;
225 if (pextent[n] <= 0)
226 {
227 pempty = 1;
228 pextent[n] = 0;
229 }
230
231 if (psize == pstride[n])
232 psize *= pextent[n];
233 else
234 psize = 0;
235 }
236 pptr = pad->data;
237 }
238 else
239 {
240 pdim = 0;
241 psize = 1;
242 pempty = 1;
243 pptr = NULL;
244 }
245
246 if (rsize != 0 && ssize != 0 && psize != 0)
247 {
248 rsize *= sizeof (GFC_REAL_4);
249 ssize *= sizeof (GFC_REAL_4);
250 psize *= sizeof (GFC_REAL_4);
251 reshape_packed ((char *)ret->data, rsize, (char *)source->data,
252 ssize, pad ? (char *)pad->data : NULL, psize);
253 return;
254 }
255 rptr = ret->data;
256 src = sptr = source->data;
257 rstride0 = rstride[0];
258 sstride0 = sstride[0];
259
260 if (sempty && pempty)
261 abort ();
262
263 if (sempty)
264 {
265 /* Pretend we are using the pad array the first time around, too. */
266 src = pptr;
267 sptr = pptr;
268 sdim = pdim;
269 for (dim = 0; dim < pdim; dim++)
270 {
271 scount[dim] = pcount[dim];
272 sextent[dim] = pextent[dim];
273 sstride[dim] = pstride[dim];
274 sstride0 = pstride[0];
275 }
276 }
277
278 while (rptr)
279 {
280 /* Select between the source and pad arrays. */
281 *rptr = *src;
282 /* Advance to the next element. */
283 rptr += rstride0;
284 src += sstride0;
285 rcount[0]++;
286 scount[0]++;
287
288 /* Advance to the next destination element. */
289 n = 0;
290 while (rcount[n] == rextent[n])
291 {
292 /* When we get to the end of a dimension, reset it and increment
293 the next dimension. */
294 rcount[n] = 0;
295 /* We could precalculate these products, but this is a less
296 frequently used path so probably not worth it. */
297 rptr -= rstride[n] * rextent[n];
298 n++;
299 if (n == rdim)
300 {
301 /* Break out of the loop. */
302 rptr = NULL;
303 break;
304 }
305 else
306 {
307 rcount[n]++;
308 rptr += rstride[n];
309 }
310 }
311 /* Advance to the next source element. */
312 n = 0;
313 while (scount[n] == sextent[n])
314 {
315 /* When we get to the end of a dimension, reset it and increment
316 the next dimension. */
317 scount[n] = 0;
318 /* We could precalculate these products, but this is a less
319 frequently used path so probably not worth it. */
320 src -= sstride[n] * sextent[n];
321 n++;
322 if (n == sdim)
323 {
324 if (sptr && pad)
325 {
326 /* Switch to the pad array. */
327 sptr = NULL;
328 sdim = pdim;
329 for (dim = 0; dim < pdim; dim++)
330 {
331 scount[dim] = pcount[dim];
332 sextent[dim] = pextent[dim];
333 sstride[dim] = pstride[dim];
334 sstride0 = sstride[0];
335 }
336 }
337 /* We now start again from the beginning of the pad array. */
338 src = pptr;
339 break;
340 }
341 else
342 {
343 scount[n]++;
344 src += sstride[n];
345 }
346 }
347 }
348 }
349
350 #endif