1 /* Single-image implementation of GNU Fortran Coarray Library
2 Copyright (C) 2011-2021 Free Software Foundation, Inc.
3 Contributed by Tobias Burnus <burnus@net-b.de>
4
5 This file is part of the GNU Fortran Coarray Runtime Library (libcaf).
6
7 Libcaf is free software; you can redistribute it and/or modify
8 it under the terms of the GNU General Public License as published by
9 the Free Software Foundation; either version 3, or (at your option)
10 any later version.
11
12 Libcaf 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 "libcaf.h"
27 #include <stdio.h> /* For fputs and fprintf. */
28 #include <stdlib.h> /* For exit and malloc. */
29 #include <string.h> /* For memcpy and memset. */
30 #include <stdarg.h> /* For variadic arguments. */
31 #include <stdint.h>
32 #include <assert.h>
33
34 /* Define GFC_CAF_CHECK to enable run-time checking. */
35 /* #define GFC_CAF_CHECK 1 */
36
37 struct caf_single_token
38 {
39 /* The pointer to the memory registered. For arrays this is the data member
40 in the descriptor. For components it's the pure data pointer. */
41 void *memptr;
42 /* The descriptor when this token is associated to an allocatable array. */
43 gfc_descriptor_t *desc;
44 /* Set when the caf lib has allocated the memory in memptr and is responsible
45 for freeing it on deregister. */
46 bool owning_memory;
47 };
48 typedef struct caf_single_token *caf_single_token_t;
49
50 #define TOKEN(X) ((caf_single_token_t) (X))
51 #define MEMTOK(X) ((caf_single_token_t) (X))->memptr
52
53 /* Single-image implementation of the CAF library.
54 Note: For performance reasons -fcoarry=single should be used
55 rather than this library. */
56
57 /* Global variables. */
58 caf_static_t *caf_static_list = NULL;
59
60 /* Keep in sync with mpi.c. */
61 static void
caf_runtime_error(const char * message,...)62 caf_runtime_error (const char *message, ...)
63 {
64 va_list ap;
65 fprintf (stderr, "Fortran runtime error: ");
66 va_start (ap, message);
67 vfprintf (stderr, message, ap);
68 va_end (ap);
69 fprintf (stderr, "\n");
70
71 /* FIXME: Shutdown the Fortran RTL to flush the buffer. PR 43849. */
72 exit (EXIT_FAILURE);
73 }
74
75 /* Error handling is similar everytime. */
76 static void
caf_internal_error(const char * msg,int * stat,char * errmsg,size_t errmsg_len,...)77 caf_internal_error (const char *msg, int *stat, char *errmsg,
78 size_t errmsg_len, ...)
79 {
80 va_list args;
81 va_start (args, errmsg_len);
82 if (stat)
83 {
84 *stat = 1;
85 if (errmsg_len > 0)
86 {
87 int len = snprintf (errmsg, errmsg_len, msg, args);
88 if (len >= 0 && errmsg_len > (size_t) len)
89 memset (&errmsg[len], ' ', errmsg_len - len);
90 }
91 va_end (args);
92 return;
93 }
94 else
95 caf_runtime_error (msg, args);
96 va_end (args);
97 }
98
99
100 void
_gfortran_caf_init(int * argc,char *** argv)101 _gfortran_caf_init (int *argc __attribute__ ((unused)),
102 char ***argv __attribute__ ((unused)))
103 {
104 }
105
106
107 void
_gfortran_caf_finalize(void)108 _gfortran_caf_finalize (void)
109 {
110 while (caf_static_list != NULL)
111 {
112 caf_static_t *tmp = caf_static_list->prev;
113 free (caf_static_list->token);
114 free (caf_static_list);
115 caf_static_list = tmp;
116 }
117 }
118
119
120 int
_gfortran_caf_this_image(int distance)121 _gfortran_caf_this_image (int distance __attribute__ ((unused)))
122 {
123 return 1;
124 }
125
126
127 int
_gfortran_caf_num_images(int distance,int failed)128 _gfortran_caf_num_images (int distance __attribute__ ((unused)),
129 int failed __attribute__ ((unused)))
130 {
131 return 1;
132 }
133
134
135 void
_gfortran_caf_register(size_t size,caf_register_t type,caf_token_t * token,gfc_descriptor_t * data,int * stat,char * errmsg,size_t errmsg_len)136 _gfortran_caf_register (size_t size, caf_register_t type, caf_token_t *token,
137 gfc_descriptor_t *data, int *stat, char *errmsg,
138 size_t errmsg_len)
139 {
140 const char alloc_fail_msg[] = "Failed to allocate coarray";
141 void *local;
142 caf_single_token_t single_token;
143
144 if (type == CAF_REGTYPE_LOCK_STATIC || type == CAF_REGTYPE_LOCK_ALLOC
145 || type == CAF_REGTYPE_CRITICAL)
146 local = calloc (size, sizeof (bool));
147 else if (type == CAF_REGTYPE_EVENT_STATIC || type == CAF_REGTYPE_EVENT_ALLOC)
148 /* In the event_(wait|post) function the counter for events is a uint32,
149 so better allocate enough memory here. */
150 local = calloc (size, sizeof (uint32_t));
151 else if (type == CAF_REGTYPE_COARRAY_ALLOC_REGISTER_ONLY)
152 local = NULL;
153 else
154 local = malloc (size);
155
156 if (type != CAF_REGTYPE_COARRAY_ALLOC_ALLOCATE_ONLY)
157 *token = malloc (sizeof (struct caf_single_token));
158
159 if (unlikely (*token == NULL
160 || (local == NULL
161 && type != CAF_REGTYPE_COARRAY_ALLOC_REGISTER_ONLY)))
162 {
163 /* Freeing the memory conditionally seems pointless, but
164 caf_internal_error () may return, when a stat is given and then the
165 memory may be lost. */
166 if (local)
167 free (local);
168 if (*token)
169 free (*token);
170 caf_internal_error (alloc_fail_msg, stat, errmsg, errmsg_len);
171 return;
172 }
173
174 single_token = TOKEN (*token);
175 single_token->memptr = local;
176 single_token->owning_memory = type != CAF_REGTYPE_COARRAY_ALLOC_REGISTER_ONLY;
177 single_token->desc = GFC_DESCRIPTOR_RANK (data) > 0 ? data : NULL;
178
179
180 if (stat)
181 *stat = 0;
182
183 if (type == CAF_REGTYPE_COARRAY_STATIC || type == CAF_REGTYPE_LOCK_STATIC
184 || type == CAF_REGTYPE_CRITICAL || type == CAF_REGTYPE_EVENT_STATIC
185 || type == CAF_REGTYPE_EVENT_ALLOC)
186 {
187 caf_static_t *tmp = malloc (sizeof (caf_static_t));
188 tmp->prev = caf_static_list;
189 tmp->token = *token;
190 caf_static_list = tmp;
191 }
192 GFC_DESCRIPTOR_DATA (data) = local;
193 }
194
195
196 void
_gfortran_caf_deregister(caf_token_t * token,caf_deregister_t type,int * stat,char * errmsg,size_t errmsg_len)197 _gfortran_caf_deregister (caf_token_t *token, caf_deregister_t type, int *stat,
198 char *errmsg __attribute__ ((unused)),
199 size_t errmsg_len __attribute__ ((unused)))
200 {
201 caf_single_token_t single_token = TOKEN (*token);
202
203 if (single_token->owning_memory && single_token->memptr)
204 free (single_token->memptr);
205
206 if (type != CAF_DEREGTYPE_COARRAY_DEALLOCATE_ONLY)
207 {
208 free (TOKEN (*token));
209 *token = NULL;
210 }
211 else
212 {
213 single_token->memptr = NULL;
214 single_token->owning_memory = false;
215 }
216
217 if (stat)
218 *stat = 0;
219 }
220
221
222 void
_gfortran_caf_sync_all(int * stat,char * errmsg,size_t errmsg_len)223 _gfortran_caf_sync_all (int *stat,
224 char *errmsg __attribute__ ((unused)),
225 size_t errmsg_len __attribute__ ((unused)))
226 {
227 __asm__ __volatile__ ("":::"memory");
228 if (stat)
229 *stat = 0;
230 }
231
232
233 void
_gfortran_caf_sync_memory(int * stat,char * errmsg,size_t errmsg_len)234 _gfortran_caf_sync_memory (int *stat,
235 char *errmsg __attribute__ ((unused)),
236 size_t errmsg_len __attribute__ ((unused)))
237 {
238 __asm__ __volatile__ ("":::"memory");
239 if (stat)
240 *stat = 0;
241 }
242
243
244 void
_gfortran_caf_sync_images(int count,int images[],int * stat,char * errmsg,size_t errmsg_len)245 _gfortran_caf_sync_images (int count __attribute__ ((unused)),
246 int images[] __attribute__ ((unused)),
247 int *stat,
248 char *errmsg __attribute__ ((unused)),
249 size_t errmsg_len __attribute__ ((unused)))
250 {
251 #ifdef GFC_CAF_CHECK
252 int i;
253
254 for (i = 0; i < count; i++)
255 if (images[i] != 1)
256 {
257 fprintf (stderr, "COARRAY ERROR: Invalid image index %d to SYNC "
258 "IMAGES", images[i]);
259 exit (EXIT_FAILURE);
260 }
261 #endif
262
263 __asm__ __volatile__ ("":::"memory");
264 if (stat)
265 *stat = 0;
266 }
267
268
269 void
_gfortran_caf_stop_numeric(int stop_code,bool quiet)270 _gfortran_caf_stop_numeric(int stop_code, bool quiet)
271 {
272 if (!quiet)
273 fprintf (stderr, "STOP %d\n", stop_code);
274 exit (0);
275 }
276
277
278 void
_gfortran_caf_stop_str(const char * string,size_t len,bool quiet)279 _gfortran_caf_stop_str(const char *string, size_t len, bool quiet)
280 {
281 if (!quiet)
282 {
283 fputs ("STOP ", stderr);
284 while (len--)
285 fputc (*(string++), stderr);
286 fputs ("\n", stderr);
287 }
288 exit (0);
289 }
290
291
292 void
_gfortran_caf_error_stop_str(const char * string,size_t len,bool quiet)293 _gfortran_caf_error_stop_str (const char *string, size_t len, bool quiet)
294 {
295 if (!quiet)
296 {
297 fputs ("ERROR STOP ", stderr);
298 while (len--)
299 fputc (*(string++), stderr);
300 fputs ("\n", stderr);
301 }
302 exit (1);
303 }
304
305
306 /* Reported that the program terminated because of a fail image issued.
307 Because this is a single image library, nothing else than aborting the whole
308 program can be done. */
309
_gfortran_caf_fail_image(void)310 void _gfortran_caf_fail_image (void)
311 {
312 fputs ("IMAGE FAILED!\n", stderr);
313 exit (0);
314 }
315
316
317 /* Get the status of image IMAGE. Because being the single image library all
318 other images are reported to be stopped. */
319
_gfortran_caf_image_status(int image,caf_team_t * team)320 int _gfortran_caf_image_status (int image,
321 caf_team_t * team __attribute__ ((unused)))
322 {
323 if (image == 1)
324 return 0;
325 else
326 return CAF_STAT_STOPPED_IMAGE;
327 }
328
329
330 /* Single image library. There cannot be any failed images with only one
331 image. */
332
333 void
_gfortran_caf_failed_images(gfc_descriptor_t * array,caf_team_t * team,int * kind)334 _gfortran_caf_failed_images (gfc_descriptor_t *array,
335 caf_team_t * team __attribute__ ((unused)),
336 int * kind)
337 {
338 int local_kind = kind != NULL ? *kind : 4;
339
340 array->base_addr = NULL;
341 array->dtype.type = BT_INTEGER;
342 array->dtype.elem_len = local_kind;
343 /* Setting lower_bound higher then upper_bound is what the compiler does to
344 indicate an empty array. */
345 array->dim[0].lower_bound = 0;
346 array->dim[0]._ubound = -1;
347 array->dim[0]._stride = 1;
348 array->offset = 0;
349 }
350
351
352 /* With only one image available no other images can be stopped. Therefore
353 return an empty array. */
354
355 void
_gfortran_caf_stopped_images(gfc_descriptor_t * array,caf_team_t * team,int * kind)356 _gfortran_caf_stopped_images (gfc_descriptor_t *array,
357 caf_team_t * team __attribute__ ((unused)),
358 int * kind)
359 {
360 int local_kind = kind != NULL ? *kind : 4;
361
362 array->base_addr = NULL;
363 array->dtype.type = BT_INTEGER;
364 array->dtype.elem_len = local_kind;
365 /* Setting lower_bound higher then upper_bound is what the compiler does to
366 indicate an empty array. */
367 array->dim[0].lower_bound = 0;
368 array->dim[0]._ubound = -1;
369 array->dim[0]._stride = 1;
370 array->offset = 0;
371 }
372
373
374 void
_gfortran_caf_error_stop(int error,bool quiet)375 _gfortran_caf_error_stop (int error, bool quiet)
376 {
377 if (!quiet)
378 fprintf (stderr, "ERROR STOP %d\n", error);
379 exit (error);
380 }
381
382
383 void
_gfortran_caf_co_broadcast(gfc_descriptor_t * a,int source_image,int * stat,char * errmsg,size_t errmsg_len)384 _gfortran_caf_co_broadcast (gfc_descriptor_t *a __attribute__ ((unused)),
385 int source_image __attribute__ ((unused)),
386 int *stat, char *errmsg __attribute__ ((unused)),
387 size_t errmsg_len __attribute__ ((unused)))
388 {
389 if (stat)
390 *stat = 0;
391 }
392
393 void
_gfortran_caf_co_sum(gfc_descriptor_t * a,int result_image,int * stat,char * errmsg,size_t errmsg_len)394 _gfortran_caf_co_sum (gfc_descriptor_t *a __attribute__ ((unused)),
395 int result_image __attribute__ ((unused)),
396 int *stat, char *errmsg __attribute__ ((unused)),
397 size_t errmsg_len __attribute__ ((unused)))
398 {
399 if (stat)
400 *stat = 0;
401 }
402
403 void
_gfortran_caf_co_min(gfc_descriptor_t * a,int result_image,int * stat,char * errmsg,int a_len,size_t errmsg_len)404 _gfortran_caf_co_min (gfc_descriptor_t *a __attribute__ ((unused)),
405 int result_image __attribute__ ((unused)),
406 int *stat, char *errmsg __attribute__ ((unused)),
407 int a_len __attribute__ ((unused)),
408 size_t errmsg_len __attribute__ ((unused)))
409 {
410 if (stat)
411 *stat = 0;
412 }
413
414 void
_gfortran_caf_co_max(gfc_descriptor_t * a,int result_image,int * stat,char * errmsg,int a_len,size_t errmsg_len)415 _gfortran_caf_co_max (gfc_descriptor_t *a __attribute__ ((unused)),
416 int result_image __attribute__ ((unused)),
417 int *stat, char *errmsg __attribute__ ((unused)),
418 int a_len __attribute__ ((unused)),
419 size_t errmsg_len __attribute__ ((unused)))
420 {
421 if (stat)
422 *stat = 0;
423 }
424
425
426 void
_gfortran_caf_co_reduce(gfc_descriptor_t * a,void * (* opr)(void *,void *),int opr_flags,int result_image,int * stat,char * errmsg,int a_len,size_t errmsg_len)427 _gfortran_caf_co_reduce (gfc_descriptor_t *a __attribute__ ((unused)),
428 void * (*opr) (void *, void *)
429 __attribute__ ((unused)),
430 int opr_flags __attribute__ ((unused)),
431 int result_image __attribute__ ((unused)),
432 int *stat, char *errmsg __attribute__ ((unused)),
433 int a_len __attribute__ ((unused)),
434 size_t errmsg_len __attribute__ ((unused)))
435 {
436 if (stat)
437 *stat = 0;
438 }
439
440
441 static void
assign_char4_from_char1(size_t dst_size,size_t src_size,uint32_t * dst,unsigned char * src)442 assign_char4_from_char1 (size_t dst_size, size_t src_size, uint32_t *dst,
443 unsigned char *src)
444 {
445 size_t i, n;
446 n = dst_size/4 > src_size ? src_size : dst_size/4;
447 for (i = 0; i < n; ++i)
448 dst[i] = (int32_t) src[i];
449 for (; i < dst_size/4; ++i)
450 dst[i] = (int32_t) ' ';
451 }
452
453
454 static void
assign_char1_from_char4(size_t dst_size,size_t src_size,unsigned char * dst,uint32_t * src)455 assign_char1_from_char4 (size_t dst_size, size_t src_size, unsigned char *dst,
456 uint32_t *src)
457 {
458 size_t i, n;
459 n = dst_size > src_size/4 ? src_size/4 : dst_size;
460 for (i = 0; i < n; ++i)
461 dst[i] = src[i] > UINT8_MAX ? (unsigned char) '?' : (unsigned char) src[i];
462 if (dst_size > n)
463 memset (&dst[n], ' ', dst_size - n);
464 }
465
466
467 static void
convert_type(void * dst,int dst_type,int dst_kind,void * src,int src_type,int src_kind,int * stat)468 convert_type (void *dst, int dst_type, int dst_kind, void *src, int src_type,
469 int src_kind, int *stat)
470 {
471 #ifdef HAVE_GFC_INTEGER_16
472 typedef __int128 int128t;
473 #else
474 typedef int64_t int128t;
475 #endif
476
477 #if defined(GFC_REAL_16_IS_LONG_DOUBLE)
478 typedef long double real128t;
479 typedef _Complex long double complex128t;
480 #elif defined(HAVE_GFC_REAL_16)
481 typedef _Complex float __attribute__((mode(TC))) __complex128;
482 typedef __float128 real128t;
483 typedef __complex128 complex128t;
484 #elif defined(HAVE_GFC_REAL_10)
485 typedef long double real128t;
486 typedef long double complex128t;
487 #else
488 typedef double real128t;
489 typedef _Complex double complex128t;
490 #endif
491
492 int128t int_val = 0;
493 real128t real_val = 0;
494 complex128t cmpx_val = 0;
495
496 switch (src_type)
497 {
498 case BT_INTEGER:
499 if (src_kind == 1)
500 int_val = *(int8_t*) src;
501 else if (src_kind == 2)
502 int_val = *(int16_t*) src;
503 else if (src_kind == 4)
504 int_val = *(int32_t*) src;
505 else if (src_kind == 8)
506 int_val = *(int64_t*) src;
507 #ifdef HAVE_GFC_INTEGER_16
508 else if (src_kind == 16)
509 int_val = *(int128t*) src;
510 #endif
511 else
512 goto error;
513 break;
514 case BT_REAL:
515 if (src_kind == 4)
516 real_val = *(float*) src;
517 else if (src_kind == 8)
518 real_val = *(double*) src;
519 #ifdef HAVE_GFC_REAL_10
520 else if (src_kind == 10)
521 real_val = *(long double*) src;
522 #endif
523 #ifdef HAVE_GFC_REAL_16
524 else if (src_kind == 16)
525 real_val = *(real128t*) src;
526 #endif
527 else
528 goto error;
529 break;
530 case BT_COMPLEX:
531 if (src_kind == 4)
532 cmpx_val = *(_Complex float*) src;
533 else if (src_kind == 8)
534 cmpx_val = *(_Complex double*) src;
535 #ifdef HAVE_GFC_REAL_10
536 else if (src_kind == 10)
537 cmpx_val = *(_Complex long double*) src;
538 #endif
539 #ifdef HAVE_GFC_REAL_16
540 else if (src_kind == 16)
541 cmpx_val = *(complex128t*) src;
542 #endif
543 else
544 goto error;
545 break;
546 default:
547 goto error;
548 }
549
550 switch (dst_type)
551 {
552 case BT_INTEGER:
553 if (src_type == BT_INTEGER)
554 {
555 if (dst_kind == 1)
556 *(int8_t*) dst = (int8_t) int_val;
557 else if (dst_kind == 2)
558 *(int16_t*) dst = (int16_t) int_val;
559 else if (dst_kind == 4)
560 *(int32_t*) dst = (int32_t) int_val;
561 else if (dst_kind == 8)
562 *(int64_t*) dst = (int64_t) int_val;
563 #ifdef HAVE_GFC_INTEGER_16
564 else if (dst_kind == 16)
565 *(int128t*) dst = (int128t) int_val;
566 #endif
567 else
568 goto error;
569 }
570 else if (src_type == BT_REAL)
571 {
572 if (dst_kind == 1)
573 *(int8_t*) dst = (int8_t) real_val;
574 else if (dst_kind == 2)
575 *(int16_t*) dst = (int16_t) real_val;
576 else if (dst_kind == 4)
577 *(int32_t*) dst = (int32_t) real_val;
578 else if (dst_kind == 8)
579 *(int64_t*) dst = (int64_t) real_val;
580 #ifdef HAVE_GFC_INTEGER_16
581 else if (dst_kind == 16)
582 *(int128t*) dst = (int128t) real_val;
583 #endif
584 else
585 goto error;
586 }
587 else if (src_type == BT_COMPLEX)
588 {
589 if (dst_kind == 1)
590 *(int8_t*) dst = (int8_t) cmpx_val;
591 else if (dst_kind == 2)
592 *(int16_t*) dst = (int16_t) cmpx_val;
593 else if (dst_kind == 4)
594 *(int32_t*) dst = (int32_t) cmpx_val;
595 else if (dst_kind == 8)
596 *(int64_t*) dst = (int64_t) cmpx_val;
597 #ifdef HAVE_GFC_INTEGER_16
598 else if (dst_kind == 16)
599 *(int128t*) dst = (int128t) cmpx_val;
600 #endif
601 else
602 goto error;
603 }
604 else
605 goto error;
606 return;
607 case BT_REAL:
608 if (src_type == BT_INTEGER)
609 {
610 if (dst_kind == 4)
611 *(float*) dst = (float) int_val;
612 else if (dst_kind == 8)
613 *(double*) dst = (double) int_val;
614 #ifdef HAVE_GFC_REAL_10
615 else if (dst_kind == 10)
616 *(long double*) dst = (long double) int_val;
617 #endif
618 #ifdef HAVE_GFC_REAL_16
619 else if (dst_kind == 16)
620 *(real128t*) dst = (real128t) int_val;
621 #endif
622 else
623 goto error;
624 }
625 else if (src_type == BT_REAL)
626 {
627 if (dst_kind == 4)
628 *(float*) dst = (float) real_val;
629 else if (dst_kind == 8)
630 *(double*) dst = (double) real_val;
631 #ifdef HAVE_GFC_REAL_10
632 else if (dst_kind == 10)
633 *(long double*) dst = (long double) real_val;
634 #endif
635 #ifdef HAVE_GFC_REAL_16
636 else if (dst_kind == 16)
637 *(real128t*) dst = (real128t) real_val;
638 #endif
639 else
640 goto error;
641 }
642 else if (src_type == BT_COMPLEX)
643 {
644 if (dst_kind == 4)
645 *(float*) dst = (float) cmpx_val;
646 else if (dst_kind == 8)
647 *(double*) dst = (double) cmpx_val;
648 #ifdef HAVE_GFC_REAL_10
649 else if (dst_kind == 10)
650 *(long double*) dst = (long double) cmpx_val;
651 #endif
652 #ifdef HAVE_GFC_REAL_16
653 else if (dst_kind == 16)
654 *(real128t*) dst = (real128t) cmpx_val;
655 #endif
656 else
657 goto error;
658 }
659 return;
660 case BT_COMPLEX:
661 if (src_type == BT_INTEGER)
662 {
663 if (dst_kind == 4)
664 *(_Complex float*) dst = (_Complex float) int_val;
665 else if (dst_kind == 8)
666 *(_Complex double*) dst = (_Complex double) int_val;
667 #ifdef HAVE_GFC_REAL_10
668 else if (dst_kind == 10)
669 *(_Complex long double*) dst = (_Complex long double) int_val;
670 #endif
671 #ifdef HAVE_GFC_REAL_16
672 else if (dst_kind == 16)
673 *(complex128t*) dst = (complex128t) int_val;
674 #endif
675 else
676 goto error;
677 }
678 else if (src_type == BT_REAL)
679 {
680 if (dst_kind == 4)
681 *(_Complex float*) dst = (_Complex float) real_val;
682 else if (dst_kind == 8)
683 *(_Complex double*) dst = (_Complex double) real_val;
684 #ifdef HAVE_GFC_REAL_10
685 else if (dst_kind == 10)
686 *(_Complex long double*) dst = (_Complex long double) real_val;
687 #endif
688 #ifdef HAVE_GFC_REAL_16
689 else if (dst_kind == 16)
690 *(complex128t*) dst = (complex128t) real_val;
691 #endif
692 else
693 goto error;
694 }
695 else if (src_type == BT_COMPLEX)
696 {
697 if (dst_kind == 4)
698 *(_Complex float*) dst = (_Complex float) cmpx_val;
699 else if (dst_kind == 8)
700 *(_Complex double*) dst = (_Complex double) cmpx_val;
701 #ifdef HAVE_GFC_REAL_10
702 else if (dst_kind == 10)
703 *(_Complex long double*) dst = (_Complex long double) cmpx_val;
704 #endif
705 #ifdef HAVE_GFC_REAL_16
706 else if (dst_kind == 16)
707 *(complex128t*) dst = (complex128t) cmpx_val;
708 #endif
709 else
710 goto error;
711 }
712 else
713 goto error;
714 return;
715 default:
716 goto error;
717 }
718
719 error:
720 fprintf (stderr, "libcaf_single RUNTIME ERROR: Cannot convert type %d kind "
721 "%d to type %d kind %d\n", src_type, src_kind, dst_type, dst_kind);
722 if (stat)
723 *stat = 1;
724 else
725 abort ();
726 }
727
728
729 void
_gfortran_caf_get(caf_token_t token,size_t offset,int image_index,gfc_descriptor_t * src,caf_vector_t * src_vector,gfc_descriptor_t * dest,int src_kind,int dst_kind,bool may_require_tmp,int * stat)730 _gfortran_caf_get (caf_token_t token, size_t offset,
731 int image_index __attribute__ ((unused)),
732 gfc_descriptor_t *src,
733 caf_vector_t *src_vector __attribute__ ((unused)),
734 gfc_descriptor_t *dest, int src_kind, int dst_kind,
735 bool may_require_tmp, int *stat)
736 {
737 /* FIXME: Handle vector subscripts. */
738 size_t i, k, size;
739 int j;
740 int rank = GFC_DESCRIPTOR_RANK (dest);
741 size_t src_size = GFC_DESCRIPTOR_SIZE (src);
742 size_t dst_size = GFC_DESCRIPTOR_SIZE (dest);
743
744 if (stat)
745 *stat = 0;
746
747 if (rank == 0)
748 {
749 void *sr = (void *) ((char *) MEMTOK (token) + offset);
750 if (GFC_DESCRIPTOR_TYPE (dest) == GFC_DESCRIPTOR_TYPE (src)
751 && dst_kind == src_kind)
752 {
753 memmove (GFC_DESCRIPTOR_DATA (dest), sr,
754 dst_size > src_size ? src_size : dst_size);
755 if (GFC_DESCRIPTOR_TYPE (dest) == BT_CHARACTER && dst_size > src_size)
756 {
757 if (dst_kind == 1)
758 memset ((void*)(char*) GFC_DESCRIPTOR_DATA (dest) + src_size,
759 ' ', dst_size - src_size);
760 else /* dst_kind == 4. */
761 for (i = src_size/4; i < dst_size/4; i++)
762 ((int32_t*) GFC_DESCRIPTOR_DATA (dest))[i] = (int32_t) ' ';
763 }
764 }
765 else if (GFC_DESCRIPTOR_TYPE (dest) == BT_CHARACTER && dst_kind == 1)
766 assign_char1_from_char4 (dst_size, src_size, GFC_DESCRIPTOR_DATA (dest),
767 sr);
768 else if (GFC_DESCRIPTOR_TYPE (dest) == BT_CHARACTER)
769 assign_char4_from_char1 (dst_size, src_size, GFC_DESCRIPTOR_DATA (dest),
770 sr);
771 else
772 convert_type (GFC_DESCRIPTOR_DATA (dest), GFC_DESCRIPTOR_TYPE (dest),
773 dst_kind, sr, GFC_DESCRIPTOR_TYPE (src), src_kind, stat);
774 return;
775 }
776
777 size = 1;
778 for (j = 0; j < rank; j++)
779 {
780 ptrdiff_t dimextent = dest->dim[j]._ubound - dest->dim[j].lower_bound + 1;
781 if (dimextent < 0)
782 dimextent = 0;
783 size *= dimextent;
784 }
785
786 if (size == 0)
787 return;
788
789 if (may_require_tmp)
790 {
791 ptrdiff_t array_offset_sr, array_offset_dst;
792 void *tmp = malloc (size*src_size);
793
794 array_offset_dst = 0;
795 for (i = 0; i < size; i++)
796 {
797 ptrdiff_t array_offset_sr = 0;
798 ptrdiff_t stride = 1;
799 ptrdiff_t extent = 1;
800 for (j = 0; j < GFC_DESCRIPTOR_RANK (src)-1; j++)
801 {
802 array_offset_sr += ((i / (extent*stride))
803 % (src->dim[j]._ubound
804 - src->dim[j].lower_bound + 1))
805 * src->dim[j]._stride;
806 extent = (src->dim[j]._ubound - src->dim[j].lower_bound + 1);
807 stride = src->dim[j]._stride;
808 }
809 array_offset_sr += (i / extent) * src->dim[rank-1]._stride;
810 void *sr = (void *)((char *) MEMTOK (token) + offset
811 + array_offset_sr*GFC_DESCRIPTOR_SIZE (src));
812 memcpy ((void *) ((char *) tmp + array_offset_dst), sr, src_size);
813 array_offset_dst += src_size;
814 }
815
816 array_offset_sr = 0;
817 for (i = 0; i < size; i++)
818 {
819 ptrdiff_t array_offset_dst = 0;
820 ptrdiff_t stride = 1;
821 ptrdiff_t extent = 1;
822 for (j = 0; j < rank-1; j++)
823 {
824 array_offset_dst += ((i / (extent*stride))
825 % (dest->dim[j]._ubound
826 - dest->dim[j].lower_bound + 1))
827 * dest->dim[j]._stride;
828 extent = (dest->dim[j]._ubound - dest->dim[j].lower_bound + 1);
829 stride = dest->dim[j]._stride;
830 }
831 array_offset_dst += (i / extent) * dest->dim[rank-1]._stride;
832 void *dst = dest->base_addr
833 + array_offset_dst*GFC_DESCRIPTOR_SIZE (dest);
834 void *sr = tmp + array_offset_sr;
835
836 if (GFC_DESCRIPTOR_TYPE (dest) == GFC_DESCRIPTOR_TYPE (src)
837 && dst_kind == src_kind)
838 {
839 memmove (dst, sr, dst_size > src_size ? src_size : dst_size);
840 if (GFC_DESCRIPTOR_TYPE (dest) == BT_CHARACTER
841 && dst_size > src_size)
842 {
843 if (dst_kind == 1)
844 memset ((void*)(char*) dst + src_size, ' ',
845 dst_size-src_size);
846 else /* dst_kind == 4. */
847 for (k = src_size/4; k < dst_size/4; k++)
848 ((int32_t*) dst)[k] = (int32_t) ' ';
849 }
850 }
851 else if (GFC_DESCRIPTOR_TYPE (dest) == BT_CHARACTER && dst_kind == 1)
852 assign_char1_from_char4 (dst_size, src_size, dst, sr);
853 else if (GFC_DESCRIPTOR_TYPE (dest) == BT_CHARACTER)
854 assign_char4_from_char1 (dst_size, src_size, dst, sr);
855 else
856 convert_type (dst, GFC_DESCRIPTOR_TYPE (dest), dst_kind,
857 sr, GFC_DESCRIPTOR_TYPE (src), src_kind, stat);
858 array_offset_sr += src_size;
859 }
860
861 free (tmp);
862 return;
863 }
864
865 for (i = 0; i < size; i++)
866 {
867 ptrdiff_t array_offset_dst = 0;
868 ptrdiff_t stride = 1;
869 ptrdiff_t extent = 1;
870 for (j = 0; j < rank-1; j++)
871 {
872 array_offset_dst += ((i / (extent*stride))
873 % (dest->dim[j]._ubound
874 - dest->dim[j].lower_bound + 1))
875 * dest->dim[j]._stride;
876 extent = (dest->dim[j]._ubound - dest->dim[j].lower_bound + 1);
877 stride = dest->dim[j]._stride;
878 }
879 array_offset_dst += (i / extent) * dest->dim[rank-1]._stride;
880 void *dst = dest->base_addr + array_offset_dst*GFC_DESCRIPTOR_SIZE (dest);
881
882 ptrdiff_t array_offset_sr = 0;
883 stride = 1;
884 extent = 1;
885 for (j = 0; j < GFC_DESCRIPTOR_RANK (src)-1; j++)
886 {
887 array_offset_sr += ((i / (extent*stride))
888 % (src->dim[j]._ubound
889 - src->dim[j].lower_bound + 1))
890 * src->dim[j]._stride;
891 extent = (src->dim[j]._ubound - src->dim[j].lower_bound + 1);
892 stride = src->dim[j]._stride;
893 }
894 array_offset_sr += (i / extent) * src->dim[rank-1]._stride;
895 void *sr = (void *)((char *) MEMTOK (token) + offset
896 + array_offset_sr*GFC_DESCRIPTOR_SIZE (src));
897
898 if (GFC_DESCRIPTOR_TYPE (dest) == GFC_DESCRIPTOR_TYPE (src)
899 && dst_kind == src_kind)
900 {
901 memmove (dst, sr, dst_size > src_size ? src_size : dst_size);
902 if (GFC_DESCRIPTOR_TYPE (dest) == BT_CHARACTER && dst_size > src_size)
903 {
904 if (dst_kind == 1)
905 memset ((void*)(char*) dst + src_size, ' ', dst_size-src_size);
906 else /* dst_kind == 4. */
907 for (k = src_size/4; k < dst_size/4; k++)
908 ((int32_t*) dst)[k] = (int32_t) ' ';
909 }
910 }
911 else if (GFC_DESCRIPTOR_TYPE (dest) == BT_CHARACTER && dst_kind == 1)
912 assign_char1_from_char4 (dst_size, src_size, dst, sr);
913 else if (GFC_DESCRIPTOR_TYPE (dest) == BT_CHARACTER)
914 assign_char4_from_char1 (dst_size, src_size, dst, sr);
915 else
916 convert_type (dst, GFC_DESCRIPTOR_TYPE (dest), dst_kind,
917 sr, GFC_DESCRIPTOR_TYPE (src), src_kind, stat);
918 }
919 }
920
921
922 void
_gfortran_caf_send(caf_token_t token,size_t offset,int image_index,gfc_descriptor_t * dest,caf_vector_t * dst_vector,gfc_descriptor_t * src,int dst_kind,int src_kind,bool may_require_tmp,int * stat)923 _gfortran_caf_send (caf_token_t token, size_t offset,
924 int image_index __attribute__ ((unused)),
925 gfc_descriptor_t *dest,
926 caf_vector_t *dst_vector __attribute__ ((unused)),
927 gfc_descriptor_t *src, int dst_kind, int src_kind,
928 bool may_require_tmp, int *stat)
929 {
930 /* FIXME: Handle vector subscripts. */
931 size_t i, k, size;
932 int j;
933 int rank = GFC_DESCRIPTOR_RANK (dest);
934 size_t src_size = GFC_DESCRIPTOR_SIZE (src);
935 size_t dst_size = GFC_DESCRIPTOR_SIZE (dest);
936
937 if (stat)
938 *stat = 0;
939
940 if (rank == 0)
941 {
942 void *dst = (void *) ((char *) MEMTOK (token) + offset);
943 if (GFC_DESCRIPTOR_TYPE (dest) == GFC_DESCRIPTOR_TYPE (src)
944 && dst_kind == src_kind)
945 {
946 memmove (dst, GFC_DESCRIPTOR_DATA (src),
947 dst_size > src_size ? src_size : dst_size);
948 if (GFC_DESCRIPTOR_TYPE (dest) == BT_CHARACTER && dst_size > src_size)
949 {
950 if (dst_kind == 1)
951 memset ((void*)(char*) dst + src_size, ' ', dst_size-src_size);
952 else /* dst_kind == 4. */
953 for (i = src_size/4; i < dst_size/4; i++)
954 ((int32_t*) dst)[i] = (int32_t) ' ';
955 }
956 }
957 else if (GFC_DESCRIPTOR_TYPE (dest) == BT_CHARACTER && dst_kind == 1)
958 assign_char1_from_char4 (dst_size, src_size, dst,
959 GFC_DESCRIPTOR_DATA (src));
960 else if (GFC_DESCRIPTOR_TYPE (dest) == BT_CHARACTER)
961 assign_char4_from_char1 (dst_size, src_size, dst,
962 GFC_DESCRIPTOR_DATA (src));
963 else
964 convert_type (dst, GFC_DESCRIPTOR_TYPE (dest), dst_kind,
965 GFC_DESCRIPTOR_DATA (src), GFC_DESCRIPTOR_TYPE (src),
966 src_kind, stat);
967 return;
968 }
969
970 size = 1;
971 for (j = 0; j < rank; j++)
972 {
973 ptrdiff_t dimextent = dest->dim[j]._ubound - dest->dim[j].lower_bound + 1;
974 if (dimextent < 0)
975 dimextent = 0;
976 size *= dimextent;
977 }
978
979 if (size == 0)
980 return;
981
982 if (may_require_tmp)
983 {
984 ptrdiff_t array_offset_sr, array_offset_dst;
985 void *tmp;
986
987 if (GFC_DESCRIPTOR_RANK (src) == 0)
988 {
989 tmp = malloc (src_size);
990 memcpy (tmp, GFC_DESCRIPTOR_DATA (src), src_size);
991 }
992 else
993 {
994 tmp = malloc (size*src_size);
995 array_offset_dst = 0;
996 for (i = 0; i < size; i++)
997 {
998 ptrdiff_t array_offset_sr = 0;
999 ptrdiff_t stride = 1;
1000 ptrdiff_t extent = 1;
1001 for (j = 0; j < GFC_DESCRIPTOR_RANK (src)-1; j++)
1002 {
1003 array_offset_sr += ((i / (extent*stride))
1004 % (src->dim[j]._ubound
1005 - src->dim[j].lower_bound + 1))
1006 * src->dim[j]._stride;
1007 extent = (src->dim[j]._ubound - src->dim[j].lower_bound + 1);
1008 stride = src->dim[j]._stride;
1009 }
1010 array_offset_sr += (i / extent) * src->dim[rank-1]._stride;
1011 void *sr = (void *) ((char *) src->base_addr
1012 + array_offset_sr*GFC_DESCRIPTOR_SIZE (src));
1013 memcpy ((void *) ((char *) tmp + array_offset_dst), sr, src_size);
1014 array_offset_dst += src_size;
1015 }
1016 }
1017
1018 array_offset_sr = 0;
1019 for (i = 0; i < size; i++)
1020 {
1021 ptrdiff_t array_offset_dst = 0;
1022 ptrdiff_t stride = 1;
1023 ptrdiff_t extent = 1;
1024 for (j = 0; j < rank-1; j++)
1025 {
1026 array_offset_dst += ((i / (extent*stride))
1027 % (dest->dim[j]._ubound
1028 - dest->dim[j].lower_bound + 1))
1029 * dest->dim[j]._stride;
1030 extent = (dest->dim[j]._ubound - dest->dim[j].lower_bound + 1);
1031 stride = dest->dim[j]._stride;
1032 }
1033 array_offset_dst += (i / extent) * dest->dim[rank-1]._stride;
1034 void *dst = (void *)((char *) MEMTOK (token) + offset
1035 + array_offset_dst*GFC_DESCRIPTOR_SIZE (dest));
1036 void *sr = tmp + array_offset_sr;
1037 if (GFC_DESCRIPTOR_TYPE (dest) == GFC_DESCRIPTOR_TYPE (src)
1038 && dst_kind == src_kind)
1039 {
1040 memmove (dst, sr,
1041 dst_size > src_size ? src_size : dst_size);
1042 if (GFC_DESCRIPTOR_TYPE (dest) == BT_CHARACTER
1043 && dst_size > src_size)
1044 {
1045 if (dst_kind == 1)
1046 memset ((void*)(char*) dst + src_size, ' ',
1047 dst_size-src_size);
1048 else /* dst_kind == 4. */
1049 for (k = src_size/4; k < dst_size/4; k++)
1050 ((int32_t*) dst)[k] = (int32_t) ' ';
1051 }
1052 }
1053 else if (GFC_DESCRIPTOR_TYPE (dest) == BT_CHARACTER && dst_kind == 1)
1054 assign_char1_from_char4 (dst_size, src_size, dst, sr);
1055 else if (GFC_DESCRIPTOR_TYPE (dest) == BT_CHARACTER)
1056 assign_char4_from_char1 (dst_size, src_size, dst, sr);
1057 else
1058 convert_type (dst, GFC_DESCRIPTOR_TYPE (dest), dst_kind,
1059 sr, GFC_DESCRIPTOR_TYPE (src), src_kind, stat);
1060 if (GFC_DESCRIPTOR_RANK (src))
1061 array_offset_sr += src_size;
1062 }
1063 free (tmp);
1064 return;
1065 }
1066
1067 for (i = 0; i < size; i++)
1068 {
1069 ptrdiff_t array_offset_dst = 0;
1070 ptrdiff_t stride = 1;
1071 ptrdiff_t extent = 1;
1072 for (j = 0; j < rank-1; j++)
1073 {
1074 array_offset_dst += ((i / (extent*stride))
1075 % (dest->dim[j]._ubound
1076 - dest->dim[j].lower_bound + 1))
1077 * dest->dim[j]._stride;
1078 extent = (dest->dim[j]._ubound - dest->dim[j].lower_bound + 1);
1079 stride = dest->dim[j]._stride;
1080 }
1081 array_offset_dst += (i / extent) * dest->dim[rank-1]._stride;
1082 void *dst = (void *)((char *) MEMTOK (token) + offset
1083 + array_offset_dst*GFC_DESCRIPTOR_SIZE (dest));
1084 void *sr;
1085 if (GFC_DESCRIPTOR_RANK (src) != 0)
1086 {
1087 ptrdiff_t array_offset_sr = 0;
1088 stride = 1;
1089 extent = 1;
1090 for (j = 0; j < GFC_DESCRIPTOR_RANK (src)-1; j++)
1091 {
1092 array_offset_sr += ((i / (extent*stride))
1093 % (src->dim[j]._ubound
1094 - src->dim[j].lower_bound + 1))
1095 * src->dim[j]._stride;
1096 extent = (src->dim[j]._ubound - src->dim[j].lower_bound + 1);
1097 stride = src->dim[j]._stride;
1098 }
1099 array_offset_sr += (i / extent) * src->dim[rank-1]._stride;
1100 sr = (void *)((char *) src->base_addr
1101 + array_offset_sr*GFC_DESCRIPTOR_SIZE (src));
1102 }
1103 else
1104 sr = src->base_addr;
1105
1106 if (GFC_DESCRIPTOR_TYPE (dest) == GFC_DESCRIPTOR_TYPE (src)
1107 && dst_kind == src_kind)
1108 {
1109 memmove (dst, sr,
1110 dst_size > src_size ? src_size : dst_size);
1111 if (GFC_DESCRIPTOR_TYPE (dest) == BT_CHARACTER && dst_size > src_size)
1112 {
1113 if (dst_kind == 1)
1114 memset ((void*)(char*) dst + src_size, ' ', dst_size-src_size);
1115 else /* dst_kind == 4. */
1116 for (k = src_size/4; k < dst_size/4; k++)
1117 ((int32_t*) dst)[k] = (int32_t) ' ';
1118 }
1119 }
1120 else if (GFC_DESCRIPTOR_TYPE (dest) == BT_CHARACTER && dst_kind == 1)
1121 assign_char1_from_char4 (dst_size, src_size, dst, sr);
1122 else if (GFC_DESCRIPTOR_TYPE (dest) == BT_CHARACTER)
1123 assign_char4_from_char1 (dst_size, src_size, dst, sr);
1124 else
1125 convert_type (dst, GFC_DESCRIPTOR_TYPE (dest), dst_kind,
1126 sr, GFC_DESCRIPTOR_TYPE (src), src_kind, stat);
1127 }
1128 }
1129
1130
1131 void
_gfortran_caf_sendget(caf_token_t dst_token,size_t dst_offset,int dst_image_index,gfc_descriptor_t * dest,caf_vector_t * dst_vector,caf_token_t src_token,size_t src_offset,int src_image_index,gfc_descriptor_t * src,caf_vector_t * src_vector,int dst_kind,int src_kind,bool may_require_tmp)1132 _gfortran_caf_sendget (caf_token_t dst_token, size_t dst_offset,
1133 int dst_image_index, gfc_descriptor_t *dest,
1134 caf_vector_t *dst_vector, caf_token_t src_token,
1135 size_t src_offset,
1136 int src_image_index __attribute__ ((unused)),
1137 gfc_descriptor_t *src,
1138 caf_vector_t *src_vector __attribute__ ((unused)),
1139 int dst_kind, int src_kind, bool may_require_tmp)
1140 {
1141 /* FIXME: Handle vector subscript of 'src_vector'. */
1142 /* For a single image, src->base_addr should be the same as src_token + offset
1143 but to play save, we do it properly. */
1144 void *src_base = GFC_DESCRIPTOR_DATA (src);
1145 GFC_DESCRIPTOR_DATA (src) = (void *) ((char *) MEMTOK (src_token)
1146 + src_offset);
1147 _gfortran_caf_send (dst_token, dst_offset, dst_image_index, dest, dst_vector,
1148 src, dst_kind, src_kind, may_require_tmp, NULL);
1149 GFC_DESCRIPTOR_DATA (src) = src_base;
1150 }
1151
1152
1153 /* Emitted when a theorectically unreachable part is reached. */
1154 const char unreachable[] = "Fatal error: unreachable alternative found.\n";
1155
1156
1157 static void
copy_data(void * ds,void * sr,int dst_type,int src_type,int dst_kind,int src_kind,size_t dst_size,size_t src_size,size_t num,int * stat)1158 copy_data (void *ds, void *sr, int dst_type, int src_type,
1159 int dst_kind, int src_kind, size_t dst_size, size_t src_size,
1160 size_t num, int *stat)
1161 {
1162 size_t k;
1163 if (dst_type == src_type && dst_kind == src_kind)
1164 {
1165 memmove (ds, sr, (dst_size > src_size ? src_size : dst_size) * num);
1166 if ((dst_type == BT_CHARACTER || src_type == BT_CHARACTER)
1167 && dst_size > src_size)
1168 {
1169 if (dst_kind == 1)
1170 memset ((void*)(char*) ds + src_size, ' ', dst_size-src_size);
1171 else /* dst_kind == 4. */
1172 for (k = src_size/4; k < dst_size/4; k++)
1173 ((int32_t*) ds)[k] = (int32_t) ' ';
1174 }
1175 }
1176 else if (dst_type == BT_CHARACTER && dst_kind == 1)
1177 assign_char1_from_char4 (dst_size, src_size, ds, sr);
1178 else if (dst_type == BT_CHARACTER)
1179 assign_char4_from_char1 (dst_size, src_size, ds, sr);
1180 else
1181 for (k = 0; k < num; ++k)
1182 {
1183 convert_type (ds, dst_type, dst_kind, sr, src_type, src_kind, stat);
1184 ds += dst_size;
1185 sr += src_size;
1186 }
1187 }
1188
1189
1190 #define COMPUTE_NUM_ITEMS(num, stride, lb, ub) \
1191 do { \
1192 index_type abs_stride = (stride) > 0 ? (stride) : -(stride); \
1193 num = (stride) > 0 ? (ub) + 1 - (lb) : (lb) + 1 - (ub); \
1194 if (num <= 0 || abs_stride < 1) return; \
1195 num = (abs_stride > 1) ? (1 + (num - 1) / abs_stride) : num; \
1196 } while (0)
1197
1198
1199 static void
get_for_ref(caf_reference_t * ref,size_t * i,size_t * dst_index,caf_single_token_t single_token,gfc_descriptor_t * dst,gfc_descriptor_t * src,void * ds,void * sr,int dst_kind,int src_kind,size_t dst_dim,size_t src_dim,size_t num,int * stat,int src_type)1200 get_for_ref (caf_reference_t *ref, size_t *i, size_t *dst_index,
1201 caf_single_token_t single_token, gfc_descriptor_t *dst,
1202 gfc_descriptor_t *src, void *ds, void *sr,
1203 int dst_kind, int src_kind, size_t dst_dim, size_t src_dim,
1204 size_t num, int *stat, int src_type)
1205 {
1206 ptrdiff_t extent_src = 1, array_offset_src = 0, stride_src;
1207 size_t next_dst_dim;
1208
1209 if (unlikely (ref == NULL))
1210 /* May be we should issue an error here, because this case should not
1211 occur. */
1212 return;
1213
1214 if (ref->next == NULL)
1215 {
1216 size_t dst_size = GFC_DESCRIPTOR_SIZE (dst);
1217 ptrdiff_t array_offset_dst = 0;;
1218 size_t dst_rank = GFC_DESCRIPTOR_RANK (dst);
1219
1220 switch (ref->type)
1221 {
1222 case CAF_REF_COMPONENT:
1223 /* Because the token is always registered after the component, its
1224 offset is always greater zero. */
1225 if (ref->u.c.caf_token_offset > 0)
1226 /* Note, that sr is dereffed here. */
1227 copy_data (ds, *(void **)(sr + ref->u.c.offset),
1228 GFC_DESCRIPTOR_TYPE (dst), src_type,
1229 dst_kind, src_kind, dst_size, ref->item_size, 1, stat);
1230 else
1231 copy_data (ds, sr + ref->u.c.offset,
1232 GFC_DESCRIPTOR_TYPE (dst), src_type,
1233 dst_kind, src_kind, dst_size, ref->item_size, 1, stat);
1234 ++(*i);
1235 return;
1236 case CAF_REF_STATIC_ARRAY:
1237 /* Intentionally fall through. */
1238 case CAF_REF_ARRAY:
1239 if (ref->u.a.mode[src_dim] == CAF_ARR_REF_NONE)
1240 {
1241 for (size_t d = 0; d < dst_rank; ++d)
1242 array_offset_dst += dst_index[d];
1243 copy_data (ds + array_offset_dst * dst_size, sr,
1244 GFC_DESCRIPTOR_TYPE (dst), src_type,
1245 dst_kind, src_kind, dst_size, ref->item_size, num,
1246 stat);
1247 *i += num;
1248 return;
1249 }
1250 break;
1251 default:
1252 caf_runtime_error (unreachable);
1253 }
1254 }
1255
1256 switch (ref->type)
1257 {
1258 case CAF_REF_COMPONENT:
1259 if (ref->u.c.caf_token_offset > 0)
1260 {
1261 single_token = *(caf_single_token_t*)(sr + ref->u.c.caf_token_offset);
1262
1263 if (ref->next && ref->next->type == CAF_REF_ARRAY)
1264 src = single_token->desc;
1265 else
1266 src = NULL;
1267
1268 if (ref->next && ref->next->type == CAF_REF_COMPONENT)
1269 /* The currently ref'ed component was allocatabe (caf_token_offset
1270 > 0) and the next ref is a component, too, then the new sr has to
1271 be dereffed. (static arrays cannot be allocatable or they
1272 become an array with descriptor. */
1273 sr = *(void **)(sr + ref->u.c.offset);
1274 else
1275 sr += ref->u.c.offset;
1276
1277 get_for_ref (ref->next, i, dst_index, single_token, dst, src,
1278 ds, sr, dst_kind, src_kind, dst_dim, 0,
1279 1, stat, src_type);
1280 }
1281 else
1282 get_for_ref (ref->next, i, dst_index, single_token, dst,
1283 (gfc_descriptor_t *)(sr + ref->u.c.offset), ds,
1284 sr + ref->u.c.offset, dst_kind, src_kind, dst_dim, 0, 1,
1285 stat, src_type);
1286 return;
1287 case CAF_REF_ARRAY:
1288 if (ref->u.a.mode[src_dim] == CAF_ARR_REF_NONE)
1289 {
1290 get_for_ref (ref->next, i, dst_index, single_token, dst,
1291 src, ds, sr, dst_kind, src_kind,
1292 dst_dim, 0, 1, stat, src_type);
1293 return;
1294 }
1295 /* Only when on the left most index switch the data pointer to
1296 the array's data pointer. */
1297 if (src_dim == 0)
1298 sr = GFC_DESCRIPTOR_DATA (src);
1299 switch (ref->u.a.mode[src_dim])
1300 {
1301 case CAF_ARR_REF_VECTOR:
1302 extent_src = GFC_DIMENSION_EXTENT (src->dim[src_dim]);
1303 array_offset_src = 0;
1304 dst_index[dst_dim] = 0;
1305 for (size_t idx = 0; idx < ref->u.a.dim[src_dim].v.nvec;
1306 ++idx)
1307 {
1308 #define KINDCASE(kind, type) case kind: \
1309 array_offset_src = (((index_type) \
1310 ((type *)ref->u.a.dim[src_dim].v.vector)[idx]) \
1311 - GFC_DIMENSION_LBOUND (src->dim[src_dim])) \
1312 * GFC_DIMENSION_STRIDE (src->dim[src_dim]); \
1313 break
1314
1315 switch (ref->u.a.dim[src_dim].v.kind)
1316 {
1317 KINDCASE (1, GFC_INTEGER_1);
1318 KINDCASE (2, GFC_INTEGER_2);
1319 KINDCASE (4, GFC_INTEGER_4);
1320 #ifdef HAVE_GFC_INTEGER_8
1321 KINDCASE (8, GFC_INTEGER_8);
1322 #endif
1323 #ifdef HAVE_GFC_INTEGER_16
1324 KINDCASE (16, GFC_INTEGER_16);
1325 #endif
1326 default:
1327 caf_runtime_error (unreachable);
1328 return;
1329 }
1330 #undef KINDCASE
1331
1332 get_for_ref (ref, i, dst_index, single_token, dst, src,
1333 ds, sr + array_offset_src * ref->item_size,
1334 dst_kind, src_kind, dst_dim + 1, src_dim + 1,
1335 1, stat, src_type);
1336 dst_index[dst_dim]
1337 += GFC_DIMENSION_STRIDE (dst->dim[dst_dim]);
1338 }
1339 return;
1340 case CAF_ARR_REF_FULL:
1341 COMPUTE_NUM_ITEMS (extent_src,
1342 ref->u.a.dim[src_dim].s.stride,
1343 GFC_DIMENSION_LBOUND (src->dim[src_dim]),
1344 GFC_DIMENSION_UBOUND (src->dim[src_dim]));
1345 stride_src = src->dim[src_dim]._stride
1346 * ref->u.a.dim[src_dim].s.stride;
1347 array_offset_src = 0;
1348 dst_index[dst_dim] = 0;
1349 for (index_type idx = 0; idx < extent_src;
1350 ++idx, array_offset_src += stride_src)
1351 {
1352 get_for_ref (ref, i, dst_index, single_token, dst, src,
1353 ds, sr + array_offset_src * ref->item_size,
1354 dst_kind, src_kind, dst_dim + 1, src_dim + 1,
1355 1, stat, src_type);
1356 dst_index[dst_dim]
1357 += GFC_DIMENSION_STRIDE (dst->dim[dst_dim]);
1358 }
1359 return;
1360 case CAF_ARR_REF_RANGE:
1361 COMPUTE_NUM_ITEMS (extent_src,
1362 ref->u.a.dim[src_dim].s.stride,
1363 ref->u.a.dim[src_dim].s.start,
1364 ref->u.a.dim[src_dim].s.end);
1365 array_offset_src = (ref->u.a.dim[src_dim].s.start
1366 - GFC_DIMENSION_LBOUND (src->dim[src_dim]))
1367 * GFC_DIMENSION_STRIDE (src->dim[src_dim]);
1368 stride_src = GFC_DIMENSION_STRIDE (src->dim[src_dim])
1369 * ref->u.a.dim[src_dim].s.stride;
1370 dst_index[dst_dim] = 0;
1371 /* Increase the dst_dim only, when the src_extent is greater one
1372 or src and dst extent are both one. Don't increase when the scalar
1373 source is not present in the dst. */
1374 next_dst_dim = extent_src > 1
1375 || (GFC_DIMENSION_EXTENT (dst->dim[dst_dim]) == 1
1376 && extent_src == 1) ? (dst_dim + 1) : dst_dim;
1377 for (index_type idx = 0; idx < extent_src; ++idx)
1378 {
1379 get_for_ref (ref, i, dst_index, single_token, dst, src,
1380 ds, sr + array_offset_src * ref->item_size,
1381 dst_kind, src_kind, next_dst_dim, src_dim + 1,
1382 1, stat, src_type);
1383 dst_index[dst_dim]
1384 += GFC_DIMENSION_STRIDE (dst->dim[dst_dim]);
1385 array_offset_src += stride_src;
1386 }
1387 return;
1388 case CAF_ARR_REF_SINGLE:
1389 array_offset_src = (ref->u.a.dim[src_dim].s.start
1390 - src->dim[src_dim].lower_bound)
1391 * GFC_DIMENSION_STRIDE (src->dim[src_dim]);
1392 dst_index[dst_dim] = 0;
1393 get_for_ref (ref, i, dst_index, single_token, dst, src, ds,
1394 sr + array_offset_src * ref->item_size,
1395 dst_kind, src_kind, dst_dim, src_dim + 1, 1,
1396 stat, src_type);
1397 return;
1398 case CAF_ARR_REF_OPEN_END:
1399 COMPUTE_NUM_ITEMS (extent_src,
1400 ref->u.a.dim[src_dim].s.stride,
1401 ref->u.a.dim[src_dim].s.start,
1402 GFC_DIMENSION_UBOUND (src->dim[src_dim]));
1403 stride_src = GFC_DIMENSION_STRIDE (src->dim[src_dim])
1404 * ref->u.a.dim[src_dim].s.stride;
1405 array_offset_src = (ref->u.a.dim[src_dim].s.start
1406 - GFC_DIMENSION_LBOUND (src->dim[src_dim]))
1407 * GFC_DIMENSION_STRIDE (src->dim[src_dim]);
1408 dst_index[dst_dim] = 0;
1409 for (index_type idx = 0; idx < extent_src; ++idx)
1410 {
1411 get_for_ref (ref, i, dst_index, single_token, dst, src,
1412 ds, sr + array_offset_src * ref->item_size,
1413 dst_kind, src_kind, dst_dim + 1, src_dim + 1,
1414 1, stat, src_type);
1415 dst_index[dst_dim]
1416 += GFC_DIMENSION_STRIDE (dst->dim[dst_dim]);
1417 array_offset_src += stride_src;
1418 }
1419 return;
1420 case CAF_ARR_REF_OPEN_START:
1421 COMPUTE_NUM_ITEMS (extent_src,
1422 ref->u.a.dim[src_dim].s.stride,
1423 GFC_DIMENSION_LBOUND (src->dim[src_dim]),
1424 ref->u.a.dim[src_dim].s.end);
1425 stride_src = GFC_DIMENSION_STRIDE (src->dim[src_dim])
1426 * ref->u.a.dim[src_dim].s.stride;
1427 array_offset_src = 0;
1428 dst_index[dst_dim] = 0;
1429 for (index_type idx = 0; idx < extent_src; ++idx)
1430 {
1431 get_for_ref (ref, i, dst_index, single_token, dst, src,
1432 ds, sr + array_offset_src * ref->item_size,
1433 dst_kind, src_kind, dst_dim + 1, src_dim + 1,
1434 1, stat, src_type);
1435 dst_index[dst_dim]
1436 += GFC_DIMENSION_STRIDE (dst->dim[dst_dim]);
1437 array_offset_src += stride_src;
1438 }
1439 return;
1440 default:
1441 caf_runtime_error (unreachable);
1442 }
1443 return;
1444 case CAF_REF_STATIC_ARRAY:
1445 if (ref->u.a.mode[src_dim] == CAF_ARR_REF_NONE)
1446 {
1447 get_for_ref (ref->next, i, dst_index, single_token, dst,
1448 NULL, ds, sr, dst_kind, src_kind,
1449 dst_dim, 0, 1, stat, src_type);
1450 return;
1451 }
1452 switch (ref->u.a.mode[src_dim])
1453 {
1454 case CAF_ARR_REF_VECTOR:
1455 array_offset_src = 0;
1456 dst_index[dst_dim] = 0;
1457 for (size_t idx = 0; idx < ref->u.a.dim[src_dim].v.nvec;
1458 ++idx)
1459 {
1460 #define KINDCASE(kind, type) case kind: \
1461 array_offset_src = ((type *)ref->u.a.dim[src_dim].v.vector)[idx]; \
1462 break
1463
1464 switch (ref->u.a.dim[src_dim].v.kind)
1465 {
1466 KINDCASE (1, GFC_INTEGER_1);
1467 KINDCASE (2, GFC_INTEGER_2);
1468 KINDCASE (4, GFC_INTEGER_4);
1469 #ifdef HAVE_GFC_INTEGER_8
1470 KINDCASE (8, GFC_INTEGER_8);
1471 #endif
1472 #ifdef HAVE_GFC_INTEGER_16
1473 KINDCASE (16, GFC_INTEGER_16);
1474 #endif
1475 default:
1476 caf_runtime_error (unreachable);
1477 return;
1478 }
1479 #undef KINDCASE
1480
1481 get_for_ref (ref, i, dst_index, single_token, dst, NULL,
1482 ds, sr + array_offset_src * ref->item_size,
1483 dst_kind, src_kind, dst_dim + 1, src_dim + 1,
1484 1, stat, src_type);
1485 dst_index[dst_dim]
1486 += GFC_DIMENSION_STRIDE (dst->dim[dst_dim]);
1487 }
1488 return;
1489 case CAF_ARR_REF_FULL:
1490 dst_index[dst_dim] = 0;
1491 for (array_offset_src = 0 ;
1492 array_offset_src <= ref->u.a.dim[src_dim].s.end;
1493 array_offset_src += ref->u.a.dim[src_dim].s.stride)
1494 {
1495 get_for_ref (ref, i, dst_index, single_token, dst, NULL,
1496 ds, sr + array_offset_src * ref->item_size,
1497 dst_kind, src_kind, dst_dim + 1, src_dim + 1,
1498 1, stat, src_type);
1499 dst_index[dst_dim]
1500 += GFC_DIMENSION_STRIDE (dst->dim[dst_dim]);
1501 }
1502 return;
1503 case CAF_ARR_REF_RANGE:
1504 COMPUTE_NUM_ITEMS (extent_src,
1505 ref->u.a.dim[src_dim].s.stride,
1506 ref->u.a.dim[src_dim].s.start,
1507 ref->u.a.dim[src_dim].s.end);
1508 array_offset_src = ref->u.a.dim[src_dim].s.start;
1509 dst_index[dst_dim] = 0;
1510 for (index_type idx = 0; idx < extent_src; ++idx)
1511 {
1512 get_for_ref (ref, i, dst_index, single_token, dst, NULL,
1513 ds, sr + array_offset_src * ref->item_size,
1514 dst_kind, src_kind, dst_dim + 1, src_dim + 1,
1515 1, stat, src_type);
1516 dst_index[dst_dim]
1517 += GFC_DIMENSION_STRIDE (dst->dim[dst_dim]);
1518 array_offset_src += ref->u.a.dim[src_dim].s.stride;
1519 }
1520 return;
1521 case CAF_ARR_REF_SINGLE:
1522 array_offset_src = ref->u.a.dim[src_dim].s.start;
1523 get_for_ref (ref, i, dst_index, single_token, dst, NULL, ds,
1524 sr + array_offset_src * ref->item_size,
1525 dst_kind, src_kind, dst_dim, src_dim + 1, 1,
1526 stat, src_type);
1527 return;
1528 /* The OPEN_* are mapped to a RANGE and therefore cannot occur. */
1529 case CAF_ARR_REF_OPEN_END:
1530 case CAF_ARR_REF_OPEN_START:
1531 default:
1532 caf_runtime_error (unreachable);
1533 }
1534 return;
1535 default:
1536 caf_runtime_error (unreachable);
1537 }
1538 }
1539
1540
1541 void
_gfortran_caf_get_by_ref(caf_token_t token,int image_index,gfc_descriptor_t * dst,caf_reference_t * refs,int dst_kind,int src_kind,bool may_require_tmp,bool dst_reallocatable,int * stat,int src_type)1542 _gfortran_caf_get_by_ref (caf_token_t token,
1543 int image_index __attribute__ ((unused)),
1544 gfc_descriptor_t *dst, caf_reference_t *refs,
1545 int dst_kind, int src_kind,
1546 bool may_require_tmp __attribute__ ((unused)),
1547 bool dst_reallocatable, int *stat,
1548 int src_type)
1549 {
1550 const char vecrefunknownkind[] = "libcaf_single::caf_get_by_ref(): "
1551 "unknown kind in vector-ref.\n";
1552 const char unknownreftype[] = "libcaf_single::caf_get_by_ref(): "
1553 "unknown reference type.\n";
1554 const char unknownarrreftype[] = "libcaf_single::caf_get_by_ref(): "
1555 "unknown array reference type.\n";
1556 const char rankoutofrange[] = "libcaf_single::caf_get_by_ref(): "
1557 "rank out of range.\n";
1558 const char extentoutofrange[] = "libcaf_single::caf_get_by_ref(): "
1559 "extent out of range.\n";
1560 const char cannotallocdst[] = "libcaf_single::caf_get_by_ref(): "
1561 "cannot allocate memory.\n";
1562 const char nonallocextentmismatch[] = "libcaf_single::caf_get_by_ref(): "
1563 "extent of non-allocatable arrays mismatch (%lu != %lu).\n";
1564 const char doublearrayref[] = "libcaf_single::caf_get_by_ref(): "
1565 "two or more array part references are not supported.\n";
1566 size_t size, i;
1567 size_t dst_index[GFC_MAX_DIMENSIONS];
1568 int dst_rank = GFC_DESCRIPTOR_RANK (dst);
1569 int dst_cur_dim = 0;
1570 size_t src_size = 0;
1571 caf_single_token_t single_token = TOKEN (token);
1572 void *memptr = single_token->memptr;
1573 gfc_descriptor_t *src = single_token->desc;
1574 caf_reference_t *riter = refs;
1575 long delta;
1576 /* Reallocation of dst.data is needed (e.g., array to small). */
1577 bool realloc_needed;
1578 /* Reallocation of dst.data is required, because data is not alloced at
1579 all. */
1580 bool realloc_required;
1581 bool extent_mismatch = false;
1582 /* Set when the first non-scalar array reference is encountered. */
1583 bool in_array_ref = false;
1584 bool array_extent_fixed = false;
1585 realloc_needed = realloc_required = GFC_DESCRIPTOR_DATA (dst) == NULL;
1586
1587 assert (!realloc_needed || dst_reallocatable);
1588
1589 if (stat)
1590 *stat = 0;
1591
1592 /* Compute the size of the result. In the beginning size just counts the
1593 number of elements. */
1594 size = 1;
1595 while (riter)
1596 {
1597 switch (riter->type)
1598 {
1599 case CAF_REF_COMPONENT:
1600 if (riter->u.c.caf_token_offset)
1601 {
1602 single_token = *(caf_single_token_t*)
1603 (memptr + riter->u.c.caf_token_offset);
1604 memptr = single_token->memptr;
1605 src = single_token->desc;
1606 }
1607 else
1608 {
1609 memptr += riter->u.c.offset;
1610 /* When the next ref is an array ref, assume there is an
1611 array descriptor at memptr. Note, static arrays do not have
1612 a descriptor. */
1613 if (riter->next && riter->next->type == CAF_REF_ARRAY)
1614 src = (gfc_descriptor_t *)memptr;
1615 else
1616 src = NULL;
1617 }
1618 break;
1619 case CAF_REF_ARRAY:
1620 for (i = 0; riter->u.a.mode[i] != CAF_ARR_REF_NONE; ++i)
1621 {
1622 switch (riter->u.a.mode[i])
1623 {
1624 case CAF_ARR_REF_VECTOR:
1625 delta = riter->u.a.dim[i].v.nvec;
1626 #define KINDCASE(kind, type) case kind: \
1627 memptr += (((index_type) \
1628 ((type *)riter->u.a.dim[i].v.vector)[0]) \
1629 - GFC_DIMENSION_LBOUND (src->dim[i])) \
1630 * GFC_DIMENSION_STRIDE (src->dim[i]) \
1631 * riter->item_size; \
1632 break
1633
1634 switch (riter->u.a.dim[i].v.kind)
1635 {
1636 KINDCASE (1, GFC_INTEGER_1);
1637 KINDCASE (2, GFC_INTEGER_2);
1638 KINDCASE (4, GFC_INTEGER_4);
1639 #ifdef HAVE_GFC_INTEGER_8
1640 KINDCASE (8, GFC_INTEGER_8);
1641 #endif
1642 #ifdef HAVE_GFC_INTEGER_16
1643 KINDCASE (16, GFC_INTEGER_16);
1644 #endif
1645 default:
1646 caf_internal_error (vecrefunknownkind, stat, NULL, 0);
1647 return;
1648 }
1649 #undef KINDCASE
1650 break;
1651 case CAF_ARR_REF_FULL:
1652 COMPUTE_NUM_ITEMS (delta,
1653 riter->u.a.dim[i].s.stride,
1654 GFC_DIMENSION_LBOUND (src->dim[i]),
1655 GFC_DIMENSION_UBOUND (src->dim[i]));
1656 /* The memptr stays unchanged when ref'ing the first element
1657 in a dimension. */
1658 break;
1659 case CAF_ARR_REF_RANGE:
1660 COMPUTE_NUM_ITEMS (delta,
1661 riter->u.a.dim[i].s.stride,
1662 riter->u.a.dim[i].s.start,
1663 riter->u.a.dim[i].s.end);
1664 memptr += (riter->u.a.dim[i].s.start
1665 - GFC_DIMENSION_LBOUND (src->dim[i]))
1666 * GFC_DIMENSION_STRIDE (src->dim[i])
1667 * riter->item_size;
1668 break;
1669 case CAF_ARR_REF_SINGLE:
1670 delta = 1;
1671 memptr += (riter->u.a.dim[i].s.start
1672 - GFC_DIMENSION_LBOUND (src->dim[i]))
1673 * GFC_DIMENSION_STRIDE (src->dim[i])
1674 * riter->item_size;
1675 break;
1676 case CAF_ARR_REF_OPEN_END:
1677 COMPUTE_NUM_ITEMS (delta,
1678 riter->u.a.dim[i].s.stride,
1679 riter->u.a.dim[i].s.start,
1680 GFC_DIMENSION_UBOUND (src->dim[i]));
1681 memptr += (riter->u.a.dim[i].s.start
1682 - GFC_DIMENSION_LBOUND (src->dim[i]))
1683 * GFC_DIMENSION_STRIDE (src->dim[i])
1684 * riter->item_size;
1685 break;
1686 case CAF_ARR_REF_OPEN_START:
1687 COMPUTE_NUM_ITEMS (delta,
1688 riter->u.a.dim[i].s.stride,
1689 GFC_DIMENSION_LBOUND (src->dim[i]),
1690 riter->u.a.dim[i].s.end);
1691 /* The memptr stays unchanged when ref'ing the first element
1692 in a dimension. */
1693 break;
1694 default:
1695 caf_internal_error (unknownarrreftype, stat, NULL, 0);
1696 return;
1697 }
1698 if (delta <= 0)
1699 return;
1700 /* Check the various properties of the destination array.
1701 Is an array expected and present? */
1702 if (delta > 1 && dst_rank == 0)
1703 {
1704 /* No, an array is required, but not provided. */
1705 caf_internal_error (extentoutofrange, stat, NULL, 0);
1706 return;
1707 }
1708 /* Special mode when called by __caf_sendget_by_ref (). */
1709 if (dst_rank == -1 && GFC_DESCRIPTOR_DATA (dst) == NULL)
1710 {
1711 dst_rank = dst_cur_dim + 1;
1712 GFC_DESCRIPTOR_RANK (dst) = dst_rank;
1713 GFC_DESCRIPTOR_SIZE (dst) = dst_kind;
1714 }
1715 /* When dst is an array. */
1716 if (dst_rank > 0)
1717 {
1718 /* Check that dst_cur_dim is valid for dst. Can be
1719 superceeded only by scalar data. */
1720 if (dst_cur_dim >= dst_rank && delta != 1)
1721 {
1722 caf_internal_error (rankoutofrange, stat, NULL, 0);
1723 return;
1724 }
1725 /* Do further checks, when the source is not scalar. */
1726 else if (delta != 1)
1727 {
1728 /* Check that the extent is not scalar and we are not in
1729 an array ref for the dst side. */
1730 if (!in_array_ref)
1731 {
1732 /* Check that this is the non-scalar extent. */
1733 if (!array_extent_fixed)
1734 {
1735 /* In an array extent now. */
1736 in_array_ref = true;
1737 /* Check that we haven't skipped any scalar
1738 dimensions yet and that the dst is
1739 compatible. */
1740 if (i > 0
1741 && dst_rank == GFC_DESCRIPTOR_RANK (src))
1742 {
1743 if (dst_reallocatable)
1744 {
1745 /* Dst is reallocatable, which means that
1746 the bounds are not set. Set them. */
1747 for (dst_cur_dim= 0; dst_cur_dim < (int)i;
1748 ++dst_cur_dim)
1749 GFC_DIMENSION_SET (dst->dim[dst_cur_dim],
1750 1, 1, 1);
1751 }
1752 else
1753 dst_cur_dim = i;
1754 }
1755 /* Else press thumbs, that there are enough
1756 dimensional refs to come. Checked below. */
1757 }
1758 else
1759 {
1760 caf_internal_error (doublearrayref, stat, NULL,
1761 0);
1762 return;
1763 }
1764 }
1765 /* When the realloc is required, then no extent may have
1766 been set. */
1767 extent_mismatch = realloc_required
1768 || GFC_DESCRIPTOR_EXTENT (dst, dst_cur_dim) != delta;
1769 /* When it already known, that a realloc is needed or
1770 the extent does not match the needed one. */
1771 if (realloc_required || realloc_needed
1772 || extent_mismatch)
1773 {
1774 /* Check whether dst is reallocatable. */
1775 if (unlikely (!dst_reallocatable))
1776 {
1777 caf_internal_error (nonallocextentmismatch, stat,
1778 NULL, 0, delta,
1779 GFC_DESCRIPTOR_EXTENT (dst,
1780 dst_cur_dim));
1781 return;
1782 }
1783 /* Only report an error, when the extent needs to be
1784 modified, which is not allowed. */
1785 else if (!dst_reallocatable && extent_mismatch)
1786 {
1787 caf_internal_error (extentoutofrange, stat, NULL,
1788 0);
1789 return;
1790 }
1791 realloc_needed = true;
1792 }
1793 /* Only change the extent when it does not match. This is
1794 to prevent resetting given array bounds. */
1795 if (extent_mismatch)
1796 GFC_DIMENSION_SET (dst->dim[dst_cur_dim], 1, delta,
1797 size);
1798 }
1799
1800 /* Only increase the dim counter, when in an array ref. */
1801 if (in_array_ref && dst_cur_dim < dst_rank)
1802 ++dst_cur_dim;
1803 }
1804 size *= (index_type)delta;
1805 }
1806 if (in_array_ref)
1807 {
1808 array_extent_fixed = true;
1809 in_array_ref = false;
1810 /* Check, if we got less dimensional refs than the rank of dst
1811 expects. */
1812 assert (dst_cur_dim == GFC_DESCRIPTOR_RANK (dst));
1813 }
1814 break;
1815 case CAF_REF_STATIC_ARRAY:
1816 for (i = 0; riter->u.a.mode[i] != CAF_ARR_REF_NONE; ++i)
1817 {
1818 switch (riter->u.a.mode[i])
1819 {
1820 case CAF_ARR_REF_VECTOR:
1821 delta = riter->u.a.dim[i].v.nvec;
1822 #define KINDCASE(kind, type) case kind: \
1823 memptr += ((type *)riter->u.a.dim[i].v.vector)[0] \
1824 * riter->item_size; \
1825 break
1826
1827 switch (riter->u.a.dim[i].v.kind)
1828 {
1829 KINDCASE (1, GFC_INTEGER_1);
1830 KINDCASE (2, GFC_INTEGER_2);
1831 KINDCASE (4, GFC_INTEGER_4);
1832 #ifdef HAVE_GFC_INTEGER_8
1833 KINDCASE (8, GFC_INTEGER_8);
1834 #endif
1835 #ifdef HAVE_GFC_INTEGER_16
1836 KINDCASE (16, GFC_INTEGER_16);
1837 #endif
1838 default:
1839 caf_internal_error (vecrefunknownkind, stat, NULL, 0);
1840 return;
1841 }
1842 #undef KINDCASE
1843 break;
1844 case CAF_ARR_REF_FULL:
1845 delta = riter->u.a.dim[i].s.end / riter->u.a.dim[i].s.stride
1846 + 1;
1847 /* The memptr stays unchanged when ref'ing the first element
1848 in a dimension. */
1849 break;
1850 case CAF_ARR_REF_RANGE:
1851 COMPUTE_NUM_ITEMS (delta,
1852 riter->u.a.dim[i].s.stride,
1853 riter->u.a.dim[i].s.start,
1854 riter->u.a.dim[i].s.end);
1855 memptr += riter->u.a.dim[i].s.start
1856 * riter->u.a.dim[i].s.stride
1857 * riter->item_size;
1858 break;
1859 case CAF_ARR_REF_SINGLE:
1860 delta = 1;
1861 memptr += riter->u.a.dim[i].s.start
1862 * riter->u.a.dim[i].s.stride
1863 * riter->item_size;
1864 break;
1865 case CAF_ARR_REF_OPEN_END:
1866 /* This and OPEN_START are mapped to a RANGE and therefore
1867 cannot occur here. */
1868 case CAF_ARR_REF_OPEN_START:
1869 default:
1870 caf_internal_error (unknownarrreftype, stat, NULL, 0);
1871 return;
1872 }
1873 if (delta <= 0)
1874 return;
1875 /* Check the various properties of the destination array.
1876 Is an array expected and present? */
1877 if (delta > 1 && dst_rank == 0)
1878 {
1879 /* No, an array is required, but not provided. */
1880 caf_internal_error (extentoutofrange, stat, NULL, 0);
1881 return;
1882 }
1883 /* Special mode when called by __caf_sendget_by_ref (). */
1884 if (dst_rank == -1 && GFC_DESCRIPTOR_DATA (dst) == NULL)
1885 {
1886 dst_rank = dst_cur_dim + 1;
1887 GFC_DESCRIPTOR_RANK (dst) = dst_rank;
1888 GFC_DESCRIPTOR_SIZE (dst) = dst_kind;
1889 }
1890 /* When dst is an array. */
1891 if (dst_rank > 0)
1892 {
1893 /* Check that dst_cur_dim is valid for dst. Can be
1894 superceeded only by scalar data. */
1895 if (dst_cur_dim >= dst_rank && delta != 1)
1896 {
1897 caf_internal_error (rankoutofrange, stat, NULL, 0);
1898 return;
1899 }
1900 /* Do further checks, when the source is not scalar. */
1901 else if (delta != 1)
1902 {
1903 /* Check that the extent is not scalar and we are not in
1904 an array ref for the dst side. */
1905 if (!in_array_ref)
1906 {
1907 /* Check that this is the non-scalar extent. */
1908 if (!array_extent_fixed)
1909 {
1910 /* In an array extent now. */
1911 in_array_ref = true;
1912 /* The dst is not reallocatable, so nothing more
1913 to do, then correct the dim counter. */
1914 dst_cur_dim = i;
1915 }
1916 else
1917 {
1918 caf_internal_error (doublearrayref, stat, NULL,
1919 0);
1920 return;
1921 }
1922 }
1923 /* When the realloc is required, then no extent may have
1924 been set. */
1925 extent_mismatch = realloc_required
1926 || GFC_DESCRIPTOR_EXTENT (dst, dst_cur_dim) != delta;
1927 /* When it is already known, that a realloc is needed or
1928 the extent does not match the needed one. */
1929 if (realloc_required || realloc_needed
1930 || extent_mismatch)
1931 {
1932 /* Check whether dst is reallocatable. */
1933 if (unlikely (!dst_reallocatable))
1934 {
1935 caf_internal_error (nonallocextentmismatch, stat,
1936 NULL, 0, delta,
1937 GFC_DESCRIPTOR_EXTENT (dst,
1938 dst_cur_dim));
1939 return;
1940 }
1941 /* Only report an error, when the extent needs to be
1942 modified, which is not allowed. */
1943 else if (!dst_reallocatable && extent_mismatch)
1944 {
1945 caf_internal_error (extentoutofrange, stat, NULL,
1946 0);
1947 return;
1948 }
1949 realloc_needed = true;
1950 }
1951 /* Only change the extent when it does not match. This is
1952 to prevent resetting given array bounds. */
1953 if (extent_mismatch)
1954 GFC_DIMENSION_SET (dst->dim[dst_cur_dim], 1, delta,
1955 size);
1956 }
1957 /* Only increase the dim counter, when in an array ref. */
1958 if (in_array_ref && dst_cur_dim < dst_rank)
1959 ++dst_cur_dim;
1960 }
1961 size *= (index_type)delta;
1962 }
1963 if (in_array_ref)
1964 {
1965 array_extent_fixed = true;
1966 in_array_ref = false;
1967 /* Check, if we got less dimensional refs than the rank of dst
1968 expects. */
1969 assert (dst_cur_dim == GFC_DESCRIPTOR_RANK (dst));
1970 }
1971 break;
1972 default:
1973 caf_internal_error (unknownreftype, stat, NULL, 0);
1974 return;
1975 }
1976 src_size = riter->item_size;
1977 riter = riter->next;
1978 }
1979 if (size == 0 || src_size == 0)
1980 return;
1981 /* Postcondition:
1982 - size contains the number of elements to store in the destination array,
1983 - src_size gives the size in bytes of each item in the destination array.
1984 */
1985
1986 if (realloc_needed)
1987 {
1988 if (!array_extent_fixed)
1989 {
1990 assert (size == 1);
1991 /* Special mode when called by __caf_sendget_by_ref (). */
1992 if (dst_rank == -1 && GFC_DESCRIPTOR_DATA (dst) == NULL)
1993 {
1994 dst_rank = dst_cur_dim + 1;
1995 GFC_DESCRIPTOR_RANK (dst) = dst_rank;
1996 GFC_DESCRIPTOR_SIZE (dst) = dst_kind;
1997 }
1998 /* This can happen only, when the result is scalar. */
1999 for (dst_cur_dim = 0; dst_cur_dim < dst_rank; ++dst_cur_dim)
2000 GFC_DIMENSION_SET (dst->dim[dst_cur_dim], 1, 1, 1);
2001 }
2002
2003 GFC_DESCRIPTOR_DATA (dst) = malloc (size * GFC_DESCRIPTOR_SIZE (dst));
2004 if (unlikely (GFC_DESCRIPTOR_DATA (dst) == NULL))
2005 {
2006 caf_internal_error (cannotallocdst, stat, NULL, 0);
2007 return;
2008 }
2009 }
2010
2011 /* Reset the token. */
2012 single_token = TOKEN (token);
2013 memptr = single_token->memptr;
2014 src = single_token->desc;
2015 memset(dst_index, 0, sizeof (dst_index));
2016 i = 0;
2017 get_for_ref (refs, &i, dst_index, single_token, dst, src,
2018 GFC_DESCRIPTOR_DATA (dst), memptr, dst_kind, src_kind, 0, 0,
2019 1, stat, src_type);
2020 }
2021
2022
2023 static void
send_by_ref(caf_reference_t * ref,size_t * i,size_t * src_index,caf_single_token_t single_token,gfc_descriptor_t * dst,gfc_descriptor_t * src,void * ds,void * sr,int dst_kind,int src_kind,size_t dst_dim,size_t src_dim,size_t num,size_t size,int * stat,int dst_type)2024 send_by_ref (caf_reference_t *ref, size_t *i, size_t *src_index,
2025 caf_single_token_t single_token, gfc_descriptor_t *dst,
2026 gfc_descriptor_t *src, void *ds, void *sr,
2027 int dst_kind, int src_kind, size_t dst_dim, size_t src_dim,
2028 size_t num, size_t size, int *stat, int dst_type)
2029 {
2030 const char vecrefunknownkind[] = "libcaf_single::caf_send_by_ref(): "
2031 "unknown kind in vector-ref.\n";
2032 ptrdiff_t extent_dst = 1, array_offset_dst = 0, stride_dst;
2033 const size_t src_rank = GFC_DESCRIPTOR_RANK (src);
2034
2035 if (unlikely (ref == NULL))
2036 /* May be we should issue an error here, because this case should not
2037 occur. */
2038 return;
2039
2040 if (ref->next == NULL)
2041 {
2042 size_t src_size = GFC_DESCRIPTOR_SIZE (src);
2043 ptrdiff_t array_offset_src = 0;;
2044
2045 switch (ref->type)
2046 {
2047 case CAF_REF_COMPONENT:
2048 if (ref->u.c.caf_token_offset > 0)
2049 {
2050 if (*(void**)(ds + ref->u.c.offset) == NULL)
2051 {
2052 /* Create a scalar temporary array descriptor. */
2053 gfc_descriptor_t static_dst;
2054 GFC_DESCRIPTOR_DATA (&static_dst) = NULL;
2055 GFC_DESCRIPTOR_DTYPE (&static_dst)
2056 = GFC_DESCRIPTOR_DTYPE (src);
2057 /* The component can be allocated now, because it is a
2058 scalar. */
2059 _gfortran_caf_register (ref->item_size,
2060 CAF_REGTYPE_COARRAY_ALLOC,
2061 ds + ref->u.c.caf_token_offset,
2062 &static_dst, stat, NULL, 0);
2063 single_token = *(caf_single_token_t *)
2064 (ds + ref->u.c.caf_token_offset);
2065 /* In case of an error in allocation return. When stat is
2066 NULL, then register_component() terminates on error. */
2067 if (stat != NULL && *stat)
2068 return;
2069 /* Publish the allocated memory. */
2070 *((void **)(ds + ref->u.c.offset))
2071 = GFC_DESCRIPTOR_DATA (&static_dst);
2072 ds = GFC_DESCRIPTOR_DATA (&static_dst);
2073 /* Set the type from the src. */
2074 dst_type = GFC_DESCRIPTOR_TYPE (src);
2075 }
2076 else
2077 {
2078 single_token = *(caf_single_token_t *)
2079 (ds + ref->u.c.caf_token_offset);
2080 dst = single_token->desc;
2081 if (dst)
2082 {
2083 ds = GFC_DESCRIPTOR_DATA (dst);
2084 dst_type = GFC_DESCRIPTOR_TYPE (dst);
2085 }
2086 else
2087 ds = *(void **)(ds + ref->u.c.offset);
2088 }
2089 copy_data (ds, sr, dst_type, GFC_DESCRIPTOR_TYPE (src),
2090 dst_kind, src_kind, ref->item_size, src_size, 1, stat);
2091 }
2092 else
2093 copy_data (ds + ref->u.c.offset, sr, dst_type,
2094 GFC_DESCRIPTOR_TYPE (src),
2095 dst_kind, src_kind, ref->item_size, src_size, 1, stat);
2096 ++(*i);
2097 return;
2098 case CAF_REF_STATIC_ARRAY:
2099 /* Intentionally fall through. */
2100 case CAF_REF_ARRAY:
2101 if (ref->u.a.mode[dst_dim] == CAF_ARR_REF_NONE)
2102 {
2103 if (src_rank > 0)
2104 {
2105 for (size_t d = 0; d < src_rank; ++d)
2106 array_offset_src += src_index[d];
2107 copy_data (ds, sr + array_offset_src * src_size,
2108 dst_type, GFC_DESCRIPTOR_TYPE (src), dst_kind,
2109 src_kind, ref->item_size, src_size, num, stat);
2110 }
2111 else
2112 copy_data (ds, sr, dst_type, GFC_DESCRIPTOR_TYPE (src),
2113 dst_kind, src_kind, ref->item_size, src_size, num,
2114 stat);
2115 *i += num;
2116 return;
2117 }
2118 break;
2119 default:
2120 caf_runtime_error (unreachable);
2121 }
2122 }
2123
2124 switch (ref->type)
2125 {
2126 case CAF_REF_COMPONENT:
2127 if (ref->u.c.caf_token_offset > 0)
2128 {
2129 if (*(void**)(ds + ref->u.c.offset) == NULL)
2130 {
2131 /* This component refs an unallocated array. Non-arrays are
2132 caught in the if (!ref->next) above. */
2133 dst = (gfc_descriptor_t *)(ds + ref->u.c.offset);
2134 /* Assume that the rank and the dimensions fit for copying src
2135 to dst. */
2136 GFC_DESCRIPTOR_DTYPE (dst) = GFC_DESCRIPTOR_DTYPE (src);
2137 dst->offset = 0;
2138 stride_dst = 1;
2139 for (size_t d = 0; d < src_rank; ++d)
2140 {
2141 extent_dst = GFC_DIMENSION_EXTENT (src->dim[d]);
2142 GFC_DIMENSION_LBOUND (dst->dim[d]) = 0;
2143 GFC_DIMENSION_UBOUND (dst->dim[d]) = extent_dst - 1;
2144 GFC_DIMENSION_STRIDE (dst->dim[d]) = stride_dst;
2145 stride_dst *= extent_dst;
2146 }
2147 /* Null the data-pointer to make register_component allocate
2148 its own memory. */
2149 GFC_DESCRIPTOR_DATA (dst) = NULL;
2150
2151 /* The size of the array is given by size. */
2152 _gfortran_caf_register (size * ref->item_size,
2153 CAF_REGTYPE_COARRAY_ALLOC,
2154 ds + ref->u.c.caf_token_offset,
2155 dst, stat, NULL, 0);
2156 /* In case of an error in allocation return. When stat is
2157 NULL, then register_component() terminates on error. */
2158 if (stat != NULL && *stat)
2159 return;
2160 }
2161 single_token = *(caf_single_token_t*)(ds + ref->u.c.caf_token_offset);
2162 /* When a component is allocatable (caf_token_offset != 0) and not an
2163 array (ref->next->type == CAF_REF_COMPONENT), then ds has to be
2164 dereffed. */
2165 if (ref->next && ref->next->type == CAF_REF_COMPONENT)
2166 ds = *(void **)(ds + ref->u.c.offset);
2167 else
2168 ds += ref->u.c.offset;
2169
2170 send_by_ref (ref->next, i, src_index, single_token,
2171 single_token->desc, src, ds, sr,
2172 dst_kind, src_kind, 0, src_dim, 1, size, stat, dst_type);
2173 }
2174 else
2175 send_by_ref (ref->next, i, src_index, single_token,
2176 (gfc_descriptor_t *)(ds + ref->u.c.offset), src,
2177 ds + ref->u.c.offset, sr, dst_kind, src_kind, 0, src_dim,
2178 1, size, stat, dst_type);
2179 return;
2180 case CAF_REF_ARRAY:
2181 if (ref->u.a.mode[dst_dim] == CAF_ARR_REF_NONE)
2182 {
2183 send_by_ref (ref->next, i, src_index, single_token,
2184 (gfc_descriptor_t *)ds, src, ds, sr, dst_kind, src_kind,
2185 0, src_dim, 1, size, stat, dst_type);
2186 return;
2187 }
2188 /* Only when on the left most index switch the data pointer to
2189 the array's data pointer. And only for non-static arrays. */
2190 if (dst_dim == 0 && ref->type != CAF_REF_STATIC_ARRAY)
2191 ds = GFC_DESCRIPTOR_DATA (dst);
2192 switch (ref->u.a.mode[dst_dim])
2193 {
2194 case CAF_ARR_REF_VECTOR:
2195 array_offset_dst = 0;
2196 src_index[src_dim] = 0;
2197 for (size_t idx = 0; idx < ref->u.a.dim[dst_dim].v.nvec;
2198 ++idx)
2199 {
2200 #define KINDCASE(kind, type) case kind: \
2201 array_offset_dst = (((index_type) \
2202 ((type *)ref->u.a.dim[dst_dim].v.vector)[idx]) \
2203 - GFC_DIMENSION_LBOUND (dst->dim[dst_dim])) \
2204 * GFC_DIMENSION_STRIDE (dst->dim[dst_dim]); \
2205 break
2206
2207 switch (ref->u.a.dim[dst_dim].v.kind)
2208 {
2209 KINDCASE (1, GFC_INTEGER_1);
2210 KINDCASE (2, GFC_INTEGER_2);
2211 KINDCASE (4, GFC_INTEGER_4);
2212 #ifdef HAVE_GFC_INTEGER_8
2213 KINDCASE (8, GFC_INTEGER_8);
2214 #endif
2215 #ifdef HAVE_GFC_INTEGER_16
2216 KINDCASE (16, GFC_INTEGER_16);
2217 #endif
2218 default:
2219 caf_internal_error (vecrefunknownkind, stat, NULL, 0);
2220 return;
2221 }
2222 #undef KINDCASE
2223
2224 send_by_ref (ref, i, src_index, single_token, dst, src,
2225 ds + array_offset_dst * ref->item_size, sr,
2226 dst_kind, src_kind, dst_dim + 1, src_dim + 1,
2227 1, size, stat, dst_type);
2228 if (src_rank > 0)
2229 src_index[src_dim]
2230 += GFC_DIMENSION_STRIDE (src->dim[src_dim]);
2231 }
2232 return;
2233 case CAF_ARR_REF_FULL:
2234 COMPUTE_NUM_ITEMS (extent_dst,
2235 ref->u.a.dim[dst_dim].s.stride,
2236 GFC_DIMENSION_LBOUND (dst->dim[dst_dim]),
2237 GFC_DIMENSION_UBOUND (dst->dim[dst_dim]));
2238 array_offset_dst = 0;
2239 stride_dst = GFC_DIMENSION_STRIDE (dst->dim[dst_dim])
2240 * ref->u.a.dim[dst_dim].s.stride;
2241 src_index[src_dim] = 0;
2242 for (index_type idx = 0; idx < extent_dst;
2243 ++idx, array_offset_dst += stride_dst)
2244 {
2245 send_by_ref (ref, i, src_index, single_token, dst, src,
2246 ds + array_offset_dst * ref->item_size, sr,
2247 dst_kind, src_kind, dst_dim + 1, src_dim + 1,
2248 1, size, stat, dst_type);
2249 if (src_rank > 0)
2250 src_index[src_dim]
2251 += GFC_DIMENSION_STRIDE (src->dim[src_dim]);
2252 }
2253 return;
2254 case CAF_ARR_REF_RANGE:
2255 COMPUTE_NUM_ITEMS (extent_dst,
2256 ref->u.a.dim[dst_dim].s.stride,
2257 ref->u.a.dim[dst_dim].s.start,
2258 ref->u.a.dim[dst_dim].s.end);
2259 array_offset_dst = ref->u.a.dim[dst_dim].s.start
2260 - GFC_DIMENSION_LBOUND (dst->dim[dst_dim]);
2261 stride_dst = GFC_DIMENSION_STRIDE (dst->dim[dst_dim])
2262 * ref->u.a.dim[dst_dim].s.stride;
2263 src_index[src_dim] = 0;
2264 for (index_type idx = 0; idx < extent_dst; ++idx)
2265 {
2266 send_by_ref (ref, i, src_index, single_token, dst, src,
2267 ds + array_offset_dst * ref->item_size, sr,
2268 dst_kind, src_kind, dst_dim + 1, src_dim + 1,
2269 1, size, stat, dst_type);
2270 if (src_rank > 0)
2271 src_index[src_dim]
2272 += GFC_DIMENSION_STRIDE (src->dim[src_dim]);
2273 array_offset_dst += stride_dst;
2274 }
2275 return;
2276 case CAF_ARR_REF_SINGLE:
2277 array_offset_dst = (ref->u.a.dim[dst_dim].s.start
2278 - GFC_DIMENSION_LBOUND (dst->dim[dst_dim]))
2279 * GFC_DIMENSION_STRIDE (dst->dim[dst_dim]);
2280 send_by_ref (ref, i, src_index, single_token, dst, src, ds
2281 + array_offset_dst * ref->item_size, sr,
2282 dst_kind, src_kind, dst_dim + 1, src_dim, 1,
2283 size, stat, dst_type);
2284 return;
2285 case CAF_ARR_REF_OPEN_END:
2286 COMPUTE_NUM_ITEMS (extent_dst,
2287 ref->u.a.dim[dst_dim].s.stride,
2288 ref->u.a.dim[dst_dim].s.start,
2289 GFC_DIMENSION_UBOUND (dst->dim[dst_dim]));
2290 array_offset_dst = ref->u.a.dim[dst_dim].s.start
2291 - GFC_DIMENSION_LBOUND (dst->dim[dst_dim]);
2292 stride_dst = GFC_DIMENSION_STRIDE (dst->dim[dst_dim])
2293 * ref->u.a.dim[dst_dim].s.stride;
2294 src_index[src_dim] = 0;
2295 for (index_type idx = 0; idx < extent_dst; ++idx)
2296 {
2297 send_by_ref (ref, i, src_index, single_token, dst, src,
2298 ds + array_offset_dst * ref->item_size, sr,
2299 dst_kind, src_kind, dst_dim + 1, src_dim + 1,
2300 1, size, stat, dst_type);
2301 if (src_rank > 0)
2302 src_index[src_dim]
2303 += GFC_DIMENSION_STRIDE (src->dim[src_dim]);
2304 array_offset_dst += stride_dst;
2305 }
2306 return;
2307 case CAF_ARR_REF_OPEN_START:
2308 COMPUTE_NUM_ITEMS (extent_dst,
2309 ref->u.a.dim[dst_dim].s.stride,
2310 GFC_DIMENSION_LBOUND (dst->dim[dst_dim]),
2311 ref->u.a.dim[dst_dim].s.end);
2312 array_offset_dst = 0;
2313 stride_dst = GFC_DIMENSION_STRIDE (dst->dim[dst_dim])
2314 * ref->u.a.dim[dst_dim].s.stride;
2315 src_index[src_dim] = 0;
2316 for (index_type idx = 0; idx < extent_dst; ++idx)
2317 {
2318 send_by_ref (ref, i, src_index, single_token, dst, src,
2319 ds + array_offset_dst * ref->item_size, sr,
2320 dst_kind, src_kind, dst_dim + 1, src_dim + 1,
2321 1, size, stat, dst_type);
2322 if (src_rank > 0)
2323 src_index[src_dim]
2324 += GFC_DIMENSION_STRIDE (src->dim[src_dim]);
2325 array_offset_dst += stride_dst;
2326 }
2327 return;
2328 default:
2329 caf_runtime_error (unreachable);
2330 }
2331 return;
2332 case CAF_REF_STATIC_ARRAY:
2333 if (ref->u.a.mode[dst_dim] == CAF_ARR_REF_NONE)
2334 {
2335 send_by_ref (ref->next, i, src_index, single_token, NULL,
2336 src, ds, sr, dst_kind, src_kind,
2337 0, src_dim, 1, size, stat, dst_type);
2338 return;
2339 }
2340 switch (ref->u.a.mode[dst_dim])
2341 {
2342 case CAF_ARR_REF_VECTOR:
2343 array_offset_dst = 0;
2344 src_index[src_dim] = 0;
2345 for (size_t idx = 0; idx < ref->u.a.dim[dst_dim].v.nvec;
2346 ++idx)
2347 {
2348 #define KINDCASE(kind, type) case kind: \
2349 array_offset_dst = ((type *)ref->u.a.dim[dst_dim].v.vector)[idx]; \
2350 break
2351
2352 switch (ref->u.a.dim[dst_dim].v.kind)
2353 {
2354 KINDCASE (1, GFC_INTEGER_1);
2355 KINDCASE (2, GFC_INTEGER_2);
2356 KINDCASE (4, GFC_INTEGER_4);
2357 #ifdef HAVE_GFC_INTEGER_8
2358 KINDCASE (8, GFC_INTEGER_8);
2359 #endif
2360 #ifdef HAVE_GFC_INTEGER_16
2361 KINDCASE (16, GFC_INTEGER_16);
2362 #endif
2363 default:
2364 caf_runtime_error (unreachable);
2365 return;
2366 }
2367 #undef KINDCASE
2368
2369 send_by_ref (ref, i, src_index, single_token, NULL, src,
2370 ds + array_offset_dst * ref->item_size, sr,
2371 dst_kind, src_kind, dst_dim + 1, src_dim + 1,
2372 1, size, stat, dst_type);
2373 src_index[src_dim]
2374 += GFC_DIMENSION_STRIDE (src->dim[src_dim]);
2375 }
2376 return;
2377 case CAF_ARR_REF_FULL:
2378 src_index[src_dim] = 0;
2379 for (array_offset_dst = 0 ;
2380 array_offset_dst <= ref->u.a.dim[dst_dim].s.end;
2381 array_offset_dst += ref->u.a.dim[dst_dim].s.stride)
2382 {
2383 send_by_ref (ref, i, src_index, single_token, NULL, src,
2384 ds + array_offset_dst * ref->item_size, sr,
2385 dst_kind, src_kind, dst_dim + 1, src_dim + 1,
2386 1, size, stat, dst_type);
2387 if (src_rank > 0)
2388 src_index[src_dim]
2389 += GFC_DIMENSION_STRIDE (src->dim[src_dim]);
2390 }
2391 return;
2392 case CAF_ARR_REF_RANGE:
2393 COMPUTE_NUM_ITEMS (extent_dst,
2394 ref->u.a.dim[dst_dim].s.stride,
2395 ref->u.a.dim[dst_dim].s.start,
2396 ref->u.a.dim[dst_dim].s.end);
2397 array_offset_dst = ref->u.a.dim[dst_dim].s.start;
2398 src_index[src_dim] = 0;
2399 for (index_type idx = 0; idx < extent_dst; ++idx)
2400 {
2401 send_by_ref (ref, i, src_index, single_token, NULL, src,
2402 ds + array_offset_dst * ref->item_size, sr,
2403 dst_kind, src_kind, dst_dim + 1, src_dim + 1,
2404 1, size, stat, dst_type);
2405 if (src_rank > 0)
2406 src_index[src_dim]
2407 += GFC_DIMENSION_STRIDE (src->dim[src_dim]);
2408 array_offset_dst += ref->u.a.dim[dst_dim].s.stride;
2409 }
2410 return;
2411 case CAF_ARR_REF_SINGLE:
2412 array_offset_dst = ref->u.a.dim[dst_dim].s.start;
2413 send_by_ref (ref, i, src_index, single_token, NULL, src,
2414 ds + array_offset_dst * ref->item_size, sr,
2415 dst_kind, src_kind, dst_dim + 1, src_dim, 1,
2416 size, stat, dst_type);
2417 return;
2418 /* The OPEN_* are mapped to a RANGE and therefore cannot occur. */
2419 case CAF_ARR_REF_OPEN_END:
2420 case CAF_ARR_REF_OPEN_START:
2421 default:
2422 caf_runtime_error (unreachable);
2423 }
2424 return;
2425 default:
2426 caf_runtime_error (unreachable);
2427 }
2428 }
2429
2430
2431 void
_gfortran_caf_send_by_ref(caf_token_t token,int image_index,gfc_descriptor_t * src,caf_reference_t * refs,int dst_kind,int src_kind,bool may_require_tmp,bool dst_reallocatable,int * stat,int dst_type)2432 _gfortran_caf_send_by_ref (caf_token_t token,
2433 int image_index __attribute__ ((unused)),
2434 gfc_descriptor_t *src, caf_reference_t *refs,
2435 int dst_kind, int src_kind,
2436 bool may_require_tmp __attribute__ ((unused)),
2437 bool dst_reallocatable, int *stat, int dst_type)
2438 {
2439 const char vecrefunknownkind[] = "libcaf_single::caf_get_by_ref(): "
2440 "unknown kind in vector-ref.\n";
2441 const char unknownreftype[] = "libcaf_single::caf_send_by_ref(): "
2442 "unknown reference type.\n";
2443 const char unknownarrreftype[] = "libcaf_single::caf_send_by_ref(): "
2444 "unknown array reference type.\n";
2445 const char rankoutofrange[] = "libcaf_single::caf_send_by_ref(): "
2446 "rank out of range.\n";
2447 const char realloconinnerref[] = "libcaf_single::caf_send_by_ref(): "
2448 "reallocation of array followed by component ref not allowed.\n";
2449 const char cannotallocdst[] = "libcaf_single::caf_send_by_ref(): "
2450 "cannot allocate memory.\n";
2451 const char nonallocextentmismatch[] = "libcaf_single::caf_send_by_ref(): "
2452 "extent of non-allocatable array mismatch.\n";
2453 const char innercompref[] = "libcaf_single::caf_send_by_ref(): "
2454 "inner unallocated component detected.\n";
2455 size_t size, i;
2456 size_t dst_index[GFC_MAX_DIMENSIONS];
2457 int src_rank = GFC_DESCRIPTOR_RANK (src);
2458 int src_cur_dim = 0;
2459 size_t src_size = 0;
2460 caf_single_token_t single_token = TOKEN (token);
2461 void *memptr = single_token->memptr;
2462 gfc_descriptor_t *dst = single_token->desc;
2463 caf_reference_t *riter = refs;
2464 long delta;
2465 bool extent_mismatch;
2466 /* Note that the component is not allocated yet. */
2467 index_type new_component_idx = -1;
2468
2469 if (stat)
2470 *stat = 0;
2471
2472 /* Compute the size of the result. In the beginning size just counts the
2473 number of elements. */
2474 size = 1;
2475 while (riter)
2476 {
2477 switch (riter->type)
2478 {
2479 case CAF_REF_COMPONENT:
2480 if (unlikely (new_component_idx != -1))
2481 {
2482 /* Allocating a component in the middle of a component ref is not
2483 support. We don't know the type to allocate. */
2484 caf_internal_error (innercompref, stat, NULL, 0);
2485 return;
2486 }
2487 if (riter->u.c.caf_token_offset > 0)
2488 {
2489 /* Check whether the allocatable component is zero, then no
2490 token is present, too. The token's pointer is not cleared
2491 when the structure is initialized. */
2492 if (*(void**)(memptr + riter->u.c.offset) == NULL)
2493 {
2494 /* This component is not yet allocated. Check that it is
2495 allocatable here. */
2496 if (!dst_reallocatable)
2497 {
2498 caf_internal_error (cannotallocdst, stat, NULL, 0);
2499 return;
2500 }
2501 single_token = NULL;
2502 memptr = NULL;
2503 dst = NULL;
2504 break;
2505 }
2506 single_token = *(caf_single_token_t*)
2507 (memptr + riter->u.c.caf_token_offset);
2508 memptr += riter->u.c.offset;
2509 dst = single_token->desc;
2510 }
2511 else
2512 {
2513 /* Regular component. */
2514 memptr += riter->u.c.offset;
2515 dst = (gfc_descriptor_t *)memptr;
2516 }
2517 break;
2518 case CAF_REF_ARRAY:
2519 if (dst != NULL)
2520 memptr = GFC_DESCRIPTOR_DATA (dst);
2521 else
2522 dst = src;
2523 /* When the dst array needs to be allocated, then look at the
2524 extent of the source array in the dimension dst_cur_dim. */
2525 for (i = 0; riter->u.a.mode[i] != CAF_ARR_REF_NONE; ++i)
2526 {
2527 switch (riter->u.a.mode[i])
2528 {
2529 case CAF_ARR_REF_VECTOR:
2530 delta = riter->u.a.dim[i].v.nvec;
2531 #define KINDCASE(kind, type) case kind: \
2532 memptr += (((index_type) \
2533 ((type *)riter->u.a.dim[i].v.vector)[0]) \
2534 - GFC_DIMENSION_LBOUND (dst->dim[i])) \
2535 * GFC_DIMENSION_STRIDE (dst->dim[i]) \
2536 * riter->item_size; \
2537 break
2538
2539 switch (riter->u.a.dim[i].v.kind)
2540 {
2541 KINDCASE (1, GFC_INTEGER_1);
2542 KINDCASE (2, GFC_INTEGER_2);
2543 KINDCASE (4, GFC_INTEGER_4);
2544 #ifdef HAVE_GFC_INTEGER_8
2545 KINDCASE (8, GFC_INTEGER_8);
2546 #endif
2547 #ifdef HAVE_GFC_INTEGER_16
2548 KINDCASE (16, GFC_INTEGER_16);
2549 #endif
2550 default:
2551 caf_internal_error (vecrefunknownkind, stat, NULL, 0);
2552 return;
2553 }
2554 #undef KINDCASE
2555 break;
2556 case CAF_ARR_REF_FULL:
2557 if (dst)
2558 COMPUTE_NUM_ITEMS (delta,
2559 riter->u.a.dim[i].s.stride,
2560 GFC_DIMENSION_LBOUND (dst->dim[i]),
2561 GFC_DIMENSION_UBOUND (dst->dim[i]));
2562 else
2563 COMPUTE_NUM_ITEMS (delta,
2564 riter->u.a.dim[i].s.stride,
2565 GFC_DIMENSION_LBOUND (src->dim[src_cur_dim]),
2566 GFC_DIMENSION_UBOUND (src->dim[src_cur_dim]));
2567 break;
2568 case CAF_ARR_REF_RANGE:
2569 COMPUTE_NUM_ITEMS (delta,
2570 riter->u.a.dim[i].s.stride,
2571 riter->u.a.dim[i].s.start,
2572 riter->u.a.dim[i].s.end);
2573 memptr += (riter->u.a.dim[i].s.start
2574 - dst->dim[i].lower_bound)
2575 * GFC_DIMENSION_STRIDE (dst->dim[i])
2576 * riter->item_size;
2577 break;
2578 case CAF_ARR_REF_SINGLE:
2579 delta = 1;
2580 memptr += (riter->u.a.dim[i].s.start
2581 - dst->dim[i].lower_bound)
2582 * GFC_DIMENSION_STRIDE (dst->dim[i])
2583 * riter->item_size;
2584 break;
2585 case CAF_ARR_REF_OPEN_END:
2586 if (dst)
2587 COMPUTE_NUM_ITEMS (delta,
2588 riter->u.a.dim[i].s.stride,
2589 riter->u.a.dim[i].s.start,
2590 GFC_DIMENSION_UBOUND (dst->dim[i]));
2591 else
2592 COMPUTE_NUM_ITEMS (delta,
2593 riter->u.a.dim[i].s.stride,
2594 riter->u.a.dim[i].s.start,
2595 GFC_DIMENSION_UBOUND (src->dim[src_cur_dim]));
2596 memptr += (riter->u.a.dim[i].s.start
2597 - dst->dim[i].lower_bound)
2598 * GFC_DIMENSION_STRIDE (dst->dim[i])
2599 * riter->item_size;
2600 break;
2601 case CAF_ARR_REF_OPEN_START:
2602 if (dst)
2603 COMPUTE_NUM_ITEMS (delta,
2604 riter->u.a.dim[i].s.stride,
2605 GFC_DIMENSION_LBOUND (dst->dim[i]),
2606 riter->u.a.dim[i].s.end);
2607 else
2608 COMPUTE_NUM_ITEMS (delta,
2609 riter->u.a.dim[i].s.stride,
2610 GFC_DIMENSION_LBOUND (src->dim[src_cur_dim]),
2611 riter->u.a.dim[i].s.end);
2612 /* The memptr stays unchanged when ref'ing the first element
2613 in a dimension. */
2614 break;
2615 default:
2616 caf_internal_error (unknownarrreftype, stat, NULL, 0);
2617 return;
2618 }
2619
2620 if (delta <= 0)
2621 return;
2622 /* Check the various properties of the source array.
2623 When src is an array. */
2624 if (delta > 1 && src_rank > 0)
2625 {
2626 /* Check that src_cur_dim is valid for src. Can be
2627 superceeded only by scalar data. */
2628 if (src_cur_dim >= src_rank)
2629 {
2630 caf_internal_error (rankoutofrange, stat, NULL, 0);
2631 return;
2632 }
2633 /* Do further checks, when the source is not scalar. */
2634 else
2635 {
2636 /* When the realloc is required, then no extent may have
2637 been set. */
2638 extent_mismatch = memptr == NULL
2639 || (dst
2640 && GFC_DESCRIPTOR_EXTENT (dst, src_cur_dim)
2641 != delta);
2642 /* When it already known, that a realloc is needed or
2643 the extent does not match the needed one. */
2644 if (extent_mismatch)
2645 {
2646 /* Check whether dst is reallocatable. */
2647 if (unlikely (!dst_reallocatable))
2648 {
2649 caf_internal_error (nonallocextentmismatch, stat,
2650 NULL, 0, delta,
2651 GFC_DESCRIPTOR_EXTENT (dst,
2652 src_cur_dim));
2653 return;
2654 }
2655 /* Report error on allocatable but missing inner
2656 ref. */
2657 else if (riter->next != NULL)
2658 {
2659 caf_internal_error (realloconinnerref, stat, NULL,
2660 0);
2661 return;
2662 }
2663 }
2664 /* Only change the extent when it does not match. This is
2665 to prevent resetting given array bounds. */
2666 if (extent_mismatch)
2667 GFC_DIMENSION_SET (dst->dim[src_cur_dim], 1, delta,
2668 size);
2669 }
2670 /* Increase the dim-counter of the src only when the extent
2671 matches. */
2672 if (src_cur_dim < src_rank
2673 && GFC_DESCRIPTOR_EXTENT (src, src_cur_dim) == delta)
2674 ++src_cur_dim;
2675 }
2676 size *= (index_type)delta;
2677 }
2678 break;
2679 case CAF_REF_STATIC_ARRAY:
2680 for (i = 0; riter->u.a.mode[i] != CAF_ARR_REF_NONE; ++i)
2681 {
2682 switch (riter->u.a.mode[i])
2683 {
2684 case CAF_ARR_REF_VECTOR:
2685 delta = riter->u.a.dim[i].v.nvec;
2686 #define KINDCASE(kind, type) case kind: \
2687 memptr += ((type *)riter->u.a.dim[i].v.vector)[0] \
2688 * riter->item_size; \
2689 break
2690
2691 switch (riter->u.a.dim[i].v.kind)
2692 {
2693 KINDCASE (1, GFC_INTEGER_1);
2694 KINDCASE (2, GFC_INTEGER_2);
2695 KINDCASE (4, GFC_INTEGER_4);
2696 #ifdef HAVE_GFC_INTEGER_8
2697 KINDCASE (8, GFC_INTEGER_8);
2698 #endif
2699 #ifdef HAVE_GFC_INTEGER_16
2700 KINDCASE (16, GFC_INTEGER_16);
2701 #endif
2702 default:
2703 caf_internal_error (vecrefunknownkind, stat, NULL, 0);
2704 return;
2705 }
2706 #undef KINDCASE
2707 break;
2708 case CAF_ARR_REF_FULL:
2709 delta = riter->u.a.dim[i].s.end / riter->u.a.dim[i].s.stride
2710 + 1;
2711 /* The memptr stays unchanged when ref'ing the first element
2712 in a dimension. */
2713 break;
2714 case CAF_ARR_REF_RANGE:
2715 COMPUTE_NUM_ITEMS (delta,
2716 riter->u.a.dim[i].s.stride,
2717 riter->u.a.dim[i].s.start,
2718 riter->u.a.dim[i].s.end);
2719 memptr += riter->u.a.dim[i].s.start
2720 * riter->u.a.dim[i].s.stride
2721 * riter->item_size;
2722 break;
2723 case CAF_ARR_REF_SINGLE:
2724 delta = 1;
2725 memptr += riter->u.a.dim[i].s.start
2726 * riter->u.a.dim[i].s.stride
2727 * riter->item_size;
2728 break;
2729 case CAF_ARR_REF_OPEN_END:
2730 /* This and OPEN_START are mapped to a RANGE and therefore
2731 cannot occur here. */
2732 case CAF_ARR_REF_OPEN_START:
2733 default:
2734 caf_internal_error (unknownarrreftype, stat, NULL, 0);
2735 return;
2736 }
2737 if (delta <= 0)
2738 return;
2739 /* Check the various properties of the source array.
2740 Only when the source array is not scalar examine its
2741 properties. */
2742 if (delta > 1 && src_rank > 0)
2743 {
2744 /* Check that src_cur_dim is valid for src. Can be
2745 superceeded only by scalar data. */
2746 if (src_cur_dim >= src_rank)
2747 {
2748 caf_internal_error (rankoutofrange, stat, NULL, 0);
2749 return;
2750 }
2751 else
2752 {
2753 /* We will not be able to realloc the dst, because that's
2754 a fixed size array. */
2755 extent_mismatch = GFC_DESCRIPTOR_EXTENT (src, src_cur_dim)
2756 != delta;
2757 /* When the extent does not match the needed one we can
2758 only stop here. */
2759 if (extent_mismatch)
2760 {
2761 caf_internal_error (nonallocextentmismatch, stat,
2762 NULL, 0, delta,
2763 GFC_DESCRIPTOR_EXTENT (src,
2764 src_cur_dim));
2765 return;
2766 }
2767 }
2768 ++src_cur_dim;
2769 }
2770 size *= (index_type)delta;
2771 }
2772 break;
2773 default:
2774 caf_internal_error (unknownreftype, stat, NULL, 0);
2775 return;
2776 }
2777 src_size = riter->item_size;
2778 riter = riter->next;
2779 }
2780 if (size == 0 || src_size == 0)
2781 return;
2782 /* Postcondition:
2783 - size contains the number of elements to store in the destination array,
2784 - src_size gives the size in bytes of each item in the destination array.
2785 */
2786
2787 /* Reset the token. */
2788 single_token = TOKEN (token);
2789 memptr = single_token->memptr;
2790 dst = single_token->desc;
2791 memset (dst_index, 0, sizeof (dst_index));
2792 i = 0;
2793 send_by_ref (refs, &i, dst_index, single_token, dst, src,
2794 memptr, GFC_DESCRIPTOR_DATA (src), dst_kind, src_kind, 0, 0,
2795 1, size, stat, dst_type);
2796 assert (i == size);
2797 }
2798
2799
2800 void
_gfortran_caf_sendget_by_ref(caf_token_t dst_token,int dst_image_index,caf_reference_t * dst_refs,caf_token_t src_token,int src_image_index,caf_reference_t * src_refs,int dst_kind,int src_kind,bool may_require_tmp,int * dst_stat,int * src_stat,int dst_type,int src_type)2801 _gfortran_caf_sendget_by_ref (caf_token_t dst_token, int dst_image_index,
2802 caf_reference_t *dst_refs, caf_token_t src_token,
2803 int src_image_index,
2804 caf_reference_t *src_refs, int dst_kind,
2805 int src_kind, bool may_require_tmp, int *dst_stat,
2806 int *src_stat, int dst_type, int src_type)
2807 {
2808 GFC_FULL_ARRAY_DESCRIPTOR (GFC_MAX_DIMENSIONS, void) temp;
2809 GFC_DESCRIPTOR_DATA (&temp) = NULL;
2810 GFC_DESCRIPTOR_RANK (&temp) = -1;
2811 GFC_DESCRIPTOR_TYPE (&temp) = dst_type;
2812
2813 _gfortran_caf_get_by_ref (src_token, src_image_index,
2814 (gfc_descriptor_t *) &temp, src_refs,
2815 dst_kind, src_kind, may_require_tmp, true,
2816 src_stat, src_type);
2817
2818 if (src_stat && *src_stat != 0)
2819 return;
2820
2821 _gfortran_caf_send_by_ref (dst_token, dst_image_index,
2822 (gfc_descriptor_t *) &temp, dst_refs,
2823 dst_kind, dst_kind, may_require_tmp, true,
2824 dst_stat, dst_type);
2825 if (GFC_DESCRIPTOR_DATA (&temp))
2826 free (GFC_DESCRIPTOR_DATA (&temp));
2827 }
2828
2829
2830 void
_gfortran_caf_atomic_define(caf_token_t token,size_t offset,int image_index,void * value,int * stat,int type,int kind)2831 _gfortran_caf_atomic_define (caf_token_t token, size_t offset,
2832 int image_index __attribute__ ((unused)),
2833 void *value, int *stat,
2834 int type __attribute__ ((unused)), int kind)
2835 {
2836 assert(kind == 4);
2837
2838 uint32_t *atom = (uint32_t *) ((char *) MEMTOK (token) + offset);
2839
2840 __atomic_store (atom, (uint32_t *) value, __ATOMIC_RELAXED);
2841
2842 if (stat)
2843 *stat = 0;
2844 }
2845
2846 void
_gfortran_caf_atomic_ref(caf_token_t token,size_t offset,int image_index,void * value,int * stat,int type,int kind)2847 _gfortran_caf_atomic_ref (caf_token_t token, size_t offset,
2848 int image_index __attribute__ ((unused)),
2849 void *value, int *stat,
2850 int type __attribute__ ((unused)), int kind)
2851 {
2852 assert(kind == 4);
2853
2854 uint32_t *atom = (uint32_t *) ((char *) MEMTOK (token) + offset);
2855
2856 __atomic_load (atom, (uint32_t *) value, __ATOMIC_RELAXED);
2857
2858 if (stat)
2859 *stat = 0;
2860 }
2861
2862
2863 void
_gfortran_caf_atomic_cas(caf_token_t token,size_t offset,int image_index,void * old,void * compare,void * new_val,int * stat,int type,int kind)2864 _gfortran_caf_atomic_cas (caf_token_t token, size_t offset,
2865 int image_index __attribute__ ((unused)),
2866 void *old, void *compare, void *new_val, int *stat,
2867 int type __attribute__ ((unused)), int kind)
2868 {
2869 assert(kind == 4);
2870
2871 uint32_t *atom = (uint32_t *) ((char *) MEMTOK (token) + offset);
2872
2873 *(uint32_t *) old = *(uint32_t *) compare;
2874 (void) __atomic_compare_exchange_n (atom, (uint32_t *) old,
2875 *(uint32_t *) new_val, false,
2876 __ATOMIC_RELAXED, __ATOMIC_RELAXED);
2877 if (stat)
2878 *stat = 0;
2879 }
2880
2881
2882 void
_gfortran_caf_atomic_op(int op,caf_token_t token,size_t offset,int image_index,void * value,void * old,int * stat,int type,int kind)2883 _gfortran_caf_atomic_op (int op, caf_token_t token, size_t offset,
2884 int image_index __attribute__ ((unused)),
2885 void *value, void *old, int *stat,
2886 int type __attribute__ ((unused)), int kind)
2887 {
2888 assert(kind == 4);
2889
2890 uint32_t res;
2891 uint32_t *atom = (uint32_t *) ((char *) MEMTOK (token) + offset);
2892
2893 switch (op)
2894 {
2895 case GFC_CAF_ATOMIC_ADD:
2896 res = __atomic_fetch_add (atom, *(uint32_t *) value, __ATOMIC_RELAXED);
2897 break;
2898 case GFC_CAF_ATOMIC_AND:
2899 res = __atomic_fetch_and (atom, *(uint32_t *) value, __ATOMIC_RELAXED);
2900 break;
2901 case GFC_CAF_ATOMIC_OR:
2902 res = __atomic_fetch_or (atom, *(uint32_t *) value, __ATOMIC_RELAXED);
2903 break;
2904 case GFC_CAF_ATOMIC_XOR:
2905 res = __atomic_fetch_xor (atom, *(uint32_t *) value, __ATOMIC_RELAXED);
2906 break;
2907 default:
2908 __builtin_unreachable();
2909 }
2910
2911 if (old)
2912 *(uint32_t *) old = res;
2913
2914 if (stat)
2915 *stat = 0;
2916 }
2917
2918 void
_gfortran_caf_event_post(caf_token_t token,size_t index,int image_index,int * stat,char * errmsg,size_t errmsg_len)2919 _gfortran_caf_event_post (caf_token_t token, size_t index,
2920 int image_index __attribute__ ((unused)),
2921 int *stat, char *errmsg __attribute__ ((unused)),
2922 size_t errmsg_len __attribute__ ((unused)))
2923 {
2924 uint32_t value = 1;
2925 uint32_t *event = (uint32_t *) ((char *) MEMTOK (token) + index
2926 * sizeof (uint32_t));
2927 __atomic_fetch_add (event, (uint32_t) value, __ATOMIC_RELAXED);
2928
2929 if(stat)
2930 *stat = 0;
2931 }
2932
2933 void
_gfortran_caf_event_wait(caf_token_t token,size_t index,int until_count,int * stat,char * errmsg,size_t errmsg_len)2934 _gfortran_caf_event_wait (caf_token_t token, size_t index,
2935 int until_count, int *stat,
2936 char *errmsg __attribute__ ((unused)),
2937 size_t errmsg_len __attribute__ ((unused)))
2938 {
2939 uint32_t *event = (uint32_t *) ((char *) MEMTOK (token) + index
2940 * sizeof (uint32_t));
2941 uint32_t value = (uint32_t)-until_count;
2942 __atomic_fetch_add (event, (uint32_t) value, __ATOMIC_RELAXED);
2943
2944 if(stat)
2945 *stat = 0;
2946 }
2947
2948 void
_gfortran_caf_event_query(caf_token_t token,size_t index,int image_index,int * count,int * stat)2949 _gfortran_caf_event_query (caf_token_t token, size_t index,
2950 int image_index __attribute__ ((unused)),
2951 int *count, int *stat)
2952 {
2953 uint32_t *event = (uint32_t *) ((char *) MEMTOK (token) + index
2954 * sizeof (uint32_t));
2955 __atomic_load (event, (uint32_t *) count, __ATOMIC_RELAXED);
2956
2957 if(stat)
2958 *stat = 0;
2959 }
2960
2961 void
_gfortran_caf_lock(caf_token_t token,size_t index,int image_index,int * acquired_lock,int * stat,char * errmsg,size_t errmsg_len)2962 _gfortran_caf_lock (caf_token_t token, size_t index,
2963 int image_index __attribute__ ((unused)),
2964 int *acquired_lock, int *stat, char *errmsg,
2965 size_t errmsg_len)
2966 {
2967 const char *msg = "Already locked";
2968 bool *lock = &((bool *) MEMTOK (token))[index];
2969
2970 if (!*lock)
2971 {
2972 *lock = true;
2973 if (acquired_lock)
2974 *acquired_lock = (int) true;
2975 if (stat)
2976 *stat = 0;
2977 return;
2978 }
2979
2980 if (acquired_lock)
2981 {
2982 *acquired_lock = (int) false;
2983 if (stat)
2984 *stat = 0;
2985 return;
2986 }
2987
2988
2989 if (stat)
2990 {
2991 *stat = 1;
2992 if (errmsg_len > 0)
2993 {
2994 size_t len = (sizeof (msg) > errmsg_len) ? errmsg_len
2995 : sizeof (msg);
2996 memcpy (errmsg, msg, len);
2997 if (errmsg_len > len)
2998 memset (&errmsg[len], ' ', errmsg_len-len);
2999 }
3000 return;
3001 }
3002 _gfortran_caf_error_stop_str (msg, strlen (msg), false);
3003 }
3004
3005
3006 void
_gfortran_caf_unlock(caf_token_t token,size_t index,int image_index,int * stat,char * errmsg,size_t errmsg_len)3007 _gfortran_caf_unlock (caf_token_t token, size_t index,
3008 int image_index __attribute__ ((unused)),
3009 int *stat, char *errmsg, size_t errmsg_len)
3010 {
3011 const char *msg = "Variable is not locked";
3012 bool *lock = &((bool *) MEMTOK (token))[index];
3013
3014 if (*lock)
3015 {
3016 *lock = false;
3017 if (stat)
3018 *stat = 0;
3019 return;
3020 }
3021
3022 if (stat)
3023 {
3024 *stat = 1;
3025 if (errmsg_len > 0)
3026 {
3027 size_t len = (sizeof (msg) > errmsg_len) ? errmsg_len
3028 : sizeof (msg);
3029 memcpy (errmsg, msg, len);
3030 if (errmsg_len > len)
3031 memset (&errmsg[len], ' ', errmsg_len-len);
3032 }
3033 return;
3034 }
3035 _gfortran_caf_error_stop_str (msg, strlen (msg), false);
3036 }
3037
3038 int
_gfortran_caf_is_present(caf_token_t token,int image_index,caf_reference_t * refs)3039 _gfortran_caf_is_present (caf_token_t token,
3040 int image_index __attribute__ ((unused)),
3041 caf_reference_t *refs)
3042 {
3043 const char arraddressingnotallowed[] = "libcaf_single::caf_is_present(): "
3044 "only scalar indexes allowed.\n";
3045 const char unknownreftype[] = "libcaf_single::caf_get_by_ref(): "
3046 "unknown reference type.\n";
3047 const char unknownarrreftype[] = "libcaf_single::caf_get_by_ref(): "
3048 "unknown array reference type.\n";
3049 size_t i;
3050 caf_single_token_t single_token = TOKEN (token);
3051 void *memptr = single_token->memptr;
3052 gfc_descriptor_t *src = single_token->desc;
3053 caf_reference_t *riter = refs;
3054
3055 while (riter)
3056 {
3057 switch (riter->type)
3058 {
3059 case CAF_REF_COMPONENT:
3060 if (riter->u.c.caf_token_offset)
3061 {
3062 single_token = *(caf_single_token_t*)
3063 (memptr + riter->u.c.caf_token_offset);
3064 memptr = single_token->memptr;
3065 src = single_token->desc;
3066 }
3067 else
3068 {
3069 memptr += riter->u.c.offset;
3070 src = (gfc_descriptor_t *)memptr;
3071 }
3072 break;
3073 case CAF_REF_ARRAY:
3074 for (i = 0; riter->u.a.mode[i] != CAF_ARR_REF_NONE; ++i)
3075 {
3076 switch (riter->u.a.mode[i])
3077 {
3078 case CAF_ARR_REF_SINGLE:
3079 memptr += (riter->u.a.dim[i].s.start
3080 - GFC_DIMENSION_LBOUND (src->dim[i]))
3081 * GFC_DIMENSION_STRIDE (src->dim[i])
3082 * riter->item_size;
3083 break;
3084 case CAF_ARR_REF_FULL:
3085 /* A full array ref is allowed on the last reference only. */
3086 if (riter->next == NULL)
3087 break;
3088 /* else fall through reporting an error. */
3089 /* FALLTHROUGH */
3090 case CAF_ARR_REF_VECTOR:
3091 case CAF_ARR_REF_RANGE:
3092 case CAF_ARR_REF_OPEN_END:
3093 case CAF_ARR_REF_OPEN_START:
3094 caf_internal_error (arraddressingnotallowed, 0, NULL, 0);
3095 return 0;
3096 default:
3097 caf_internal_error (unknownarrreftype, 0, NULL, 0);
3098 return 0;
3099 }
3100 }
3101 break;
3102 case CAF_REF_STATIC_ARRAY:
3103 for (i = 0; riter->u.a.mode[i] != CAF_ARR_REF_NONE; ++i)
3104 {
3105 switch (riter->u.a.mode[i])
3106 {
3107 case CAF_ARR_REF_SINGLE:
3108 memptr += riter->u.a.dim[i].s.start
3109 * riter->u.a.dim[i].s.stride
3110 * riter->item_size;
3111 break;
3112 case CAF_ARR_REF_FULL:
3113 /* A full array ref is allowed on the last reference only. */
3114 if (riter->next == NULL)
3115 break;
3116 /* else fall through reporting an error. */
3117 /* FALLTHROUGH */
3118 case CAF_ARR_REF_VECTOR:
3119 case CAF_ARR_REF_RANGE:
3120 case CAF_ARR_REF_OPEN_END:
3121 case CAF_ARR_REF_OPEN_START:
3122 caf_internal_error (arraddressingnotallowed, 0, NULL, 0);
3123 return 0;
3124 default:
3125 caf_internal_error (unknownarrreftype, 0, NULL, 0);
3126 return 0;
3127 }
3128 }
3129 break;
3130 default:
3131 caf_internal_error (unknownreftype, 0, NULL, 0);
3132 return 0;
3133 }
3134 riter = riter->next;
3135 }
3136 return memptr != NULL;
3137 }
3138
3139 /* Reference the libraries implementation. */
3140 extern void _gfortran_random_init (int32_t, int32_t, int32_t);
3141
_gfortran_caf_random_init(bool repeatable,bool image_distinct)3142 void _gfortran_caf_random_init (bool repeatable, bool image_distinct)
3143 {
3144 /* In a single image implementation always forward to the gfortran
3145 routine. */
3146 _gfortran_random_init (repeatable, image_distinct, 1);
3147 }
3148