YetAnotherCoupler 3.2.0
Loading...
Searching...
No Matches
couple_config.c
Go to the documentation of this file.
1// Copyright (c) 2024 The YAC Authors
2//
3// SPDX-License-Identifier: BSD-3-Clause
4
5#ifdef HAVE_CONFIG_H
6#include "config.h"
7#endif
8
9#include <stdlib.h>
10#include <stdio.h>
11#include <string.h>
12#include <math.h>
13
14#include "couple_config.h"
15#include "yac_mpi_common.h"
16#include "yac.h"
17#include "mtime_calendar.h"
18#include "utils_common.h"
19#include "config_yaml.h"
20
22 char const * name;
23 char* metadata;
24};
25
27 char const * name;
28 size_t grid_idx;
29 char const * timestep;
30 char * metadata;
33};
34
36
38 size_t num_fields;
39
40 char const * name;
41 char * metadata;
42};
43
68
76
94
96
97 struct yac_couple_config * couple_config =
98 xmalloc(1 * sizeof(*couple_config));
99
100 couple_config->start_datetime = NULL;
101 couple_config->end_datetime = NULL;
102
103 couple_config->config_filename = NULL;
105
106 couple_config->grids = NULL;
107 couple_config->num_grids = 0;
108
109 couple_config->components = NULL;
110 couple_config->num_components = 0;
111
112 couple_config->couples = NULL;
113 couple_config->num_couples = 0;
114
115 return couple_config;
116}
117
118static char const * string_dup(char const * string) {
119 return (string != NULL)?strdup(string):NULL;
120}
121
122static void yac_couple_config_field_free(void * field_){
123 struct yac_couple_config_field * field = field_;
124 free((void*)field->name);
125 free((void*)field->timestep);
126 free(field->metadata);
127}
128
129
130static void yac_couple_config_component_free(void * component_) {
131
132 struct yac_couple_config_component * component = component_;
133
134 for (size_t i = 0; i < component->num_fields; ++i)
135 yac_couple_config_field_free(component->fields + i);
136 free(component->fields);
137 free((void*)(component->name));
138 free(component->metadata);
139}
140
141static void yac_couple_config_field_couple_free(void * field_couple_) {
142
143 struct yac_couple_config_field_couple * field_couple = field_couple_;
145 free((void*)field_couple->coupling_period);
146 free((void*)field_couple->weight_file_name);
147 for (size_t i = 0; i < field_couple->num_src_mask_names; ++i)
148 free((void*)field_couple->src_mask_names[i]);
149 free(field_couple->src_mask_names);
150 free(field_couple->tgt_mask_name);
151}
152
153static void yac_couple_config_couple_free(void * couple_) {
154
155 struct yac_couple_config_couple * couple = couple_;
156 for (size_t i = 0; i < couple->num_field_couples; ++i)
158 free(couple->field_couples);
159 couple->field_couples = NULL;
160 couple->num_field_couples = 0;
161}
162
163static void yac_couple_config_grid_free(void * grid_) {
164
165 struct yac_couple_config_grid * grid = grid_;
166
167 free((void*)grid->name);
168 free(grid->metadata);
169}
170
171static int yac_couple_config_grid_compare(void const * a, void const * b) {
172 YAC_ASSERT(((struct yac_couple_config_grid const *)a)->name != NULL &&
173 ((struct yac_couple_config_grid const *)b)->name != NULL,
174 "ERROR(yac_couple_config_grid_compare): "
175 "invalid name (NULL is not allowed)");
176 return strcmp(((struct yac_couple_config_grid const *)a)->name,
177 ((struct yac_couple_config_grid const *)b)->name);
178}
179
180static int yac_couple_config_field_compare(void const * a_, void const * b_) {
181 struct yac_couple_config_field const * a = a_;
182 struct yac_couple_config_field const * b = b_;
183 int ret = (a->grid_idx > b->grid_idx) - (a->grid_idx < b->grid_idx);
184 if (ret) return ret;
186 a->name != NULL && b->name != NULL,
187 "ERROR(yac_couple_config_field_compare): "
188 "invalid name (NULL is not allowed)");
189 return strcmp(a->name, b->name);
190}
191
193 void const * a_, void const * b_) {
194 struct yac_couple_config_component const * a = a_;
195 struct yac_couple_config_component const * b = b_;
196 YAC_ASSERT(a->name != NULL && b->name != NULL,
197 "ERROR(yac_couple_config_component_compare): "
198 "invalid name (NULL is not allowed)");
199 return strcmp(a->name, b->name);
200}
201
203 void const * a_, void const * b_) {
204 struct yac_couple_config_field_couple const * a = a_;
205 struct yac_couple_config_field_couple const * b = b_;
206 int ret;
207 ret = (a->source.component_idx > b->source.component_idx) -
209 if (ret) return ret;
210 ret = (a->target.component_idx > b->target.component_idx) -
212 if (ret) return ret;
213 ret = (a->source.field_idx > b->source.field_idx) -
215 if (ret) return ret;
216 return (a->target.field_idx > b->target.field_idx) -
218}
219
221 void const * a_, void const * b_) {
222 struct yac_couple_config_couple const * a = a_;
223 struct yac_couple_config_couple const * b = b_;
224 int ret;
225 ret = (a->component_indices[0] > b->component_indices[0]) -
226 (a->component_indices[0] < b->component_indices[0]);
227 if (ret) return ret;
228 return (a->component_indices[1] > b->component_indices[1]) -
229 (a->component_indices[1] < b->component_indices[1]);
230}
231
233 char const * string_name, char ** string, MPI_Comm comm) {
234
235 int rank;
236 MPI_Comm_rank(comm, &rank);
237 struct {
238 int len;
239 int rank;
240 } data_pair;
241 size_t len = (*string != NULL)?(strlen(*string) + 1):0;
243 len <= INT_MAX,
244 "ERROR(couple_config_sync_string): \"%s\" too long", string_name);
245 data_pair.len = (int)len;
246 data_pair.rank = rank;
247
248 // determine the broadcaster
250 MPI_Allreduce(
251 MPI_IN_PLACE, &data_pair, 1, MPI_2INT, MPI_MAXLOC, comm), comm);
252 if (data_pair.len == 0) return;
253
254 // broadcast string
255 char const * string_bak = NULL;
256 if (data_pair.rank != rank) {
257 string_bak = *string;
258 *string = xmalloc((size_t)data_pair.len * sizeof(**string));
259 }
261 MPI_Bcast(*string, data_pair.len, MPI_CHAR, data_pair.rank, comm), comm);
262
263 // check for consistency
264 if (data_pair.rank != rank) {
266 (string_bak == NULL) ||
267 !strcmp(string_bak, *string),
268 "ERROR(couple_config_sync_string): inconsistent \"%s\" definition "
269 "(\"%s\" != \"%s\")", string_name, string_bak, *string);
270 free((void*)string_bak);
271 }
272}
273
274static void couple_config_sync_calendar(MPI_Comm comm) {
275
276 int calendar = (int)getCalendarType();
277 if (calendar == CALENDAR_NOT_SET) calendar = INT_MAX;
278
279 // broadcast calendar
281 MPI_Allreduce(
282 MPI_IN_PLACE, &calendar, 1, MPI_INT, MPI_MIN, comm), comm);
283
284 // if no process has defined a calendar
285 if (calendar == INT_MAX) return;
286
287 // set calendar (if local process has already defined a calendar,
288 // the definition is checked for consistency)
289 yac_cdef_calendar(calendar);
290}
291
293 enum yac_config_filetype * filetype, MPI_Comm comm) {
294
295 int filetype_ = (int)*filetype;
296
297 // distribute filetype
299 MPI_Allreduce(
300 MPI_IN_PLACE, &filetype_, 1, MPI_INT, MPI_MIN, comm), comm);
301
302 // check for consistency
304 (*filetype == YAC_CONFIG_FILETYPE_INVALID) ||
305 ((int)*filetype == filetype_),
306 "ERROR(couple_config_sync_config_filetype): "
307 "inconsistent defintion of output configuration filetype");
308
309 *filetype = (enum yac_config_filetype)filetype_;
310}
311
313 void * a_, void * b_, MPI_Comm comm) {
314
315 struct yac_couple_config_component * a = a_;
316 struct yac_couple_config_component * b = b_;
317
318 if (b) {
319 a->num_fields = b->num_fields;
320 a->fields = b->fields;
321 b->num_fields = 0;
322 b->fields = NULL;
323 a->metadata = b->metadata;
324 b->metadata = NULL;
325 }
326
328 "component metadata", (char **)&(a->metadata), comm);
329}
330
332 void * a_, void * b_, MPI_Comm comm) {
333
334 struct yac_couple_config_grid * a = a_;
335 struct yac_couple_config_grid * b = b_;
336
337 if (b!= NULL) {
338 a->metadata = b->metadata;
339 b->metadata = NULL;
340 }
341
342 couple_config_sync_string("grid metadata", (char **)&(a->metadata), comm);
343}
344
346 void * a_, void * b_, MPI_Comm comm) {
347
348 struct yac_couple_config_field * a = a_;
349 struct yac_couple_config_field * b = b_;
350
351 enum {
352 TIMESTEP_IDX,
353 METADATA_IDX,
354 FRAC_MASK_IDX,
355 COLLECTION_SIZE_IDX,
356 DATA_PAIR_COUNT
357 };
358
359 int rank;
360 MPI_Comm_rank(comm, &rank);
361 struct {
362 int len;
363 int rank;
364 } data_pairs[DATA_PAIR_COUNT];
365 size_t timestep_len =
366 ((b != NULL) && (b->timestep != NULL))?(strlen(b->timestep) + 1):0;
367 size_t metadata_len =
368 ((b != NULL) && (b->metadata != NULL))?(strlen(b->metadata) + 1):0;
370 timestep_len <= INT_MAX,
371 "ERROR(yac_couple_config_field_merge): timestep string too long");
373 metadata_len <= INT_MAX,
374 "ERROR(yac_couple_config_field_merge): metadata string too long");
376 (b == NULL) || (b->collection_size <= INT_MAX) ||
377 (b->collection_size == SIZE_MAX),
378 "ERROR(yac_couple_config_field_merge): invalid collection size \"%zu\"",
379 b->collection_size);
380 data_pairs[TIMESTEP_IDX].len = (int)timestep_len;
381 data_pairs[TIMESTEP_IDX].rank = rank;
382 data_pairs[METADATA_IDX].len = (int)metadata_len;
383 data_pairs[METADATA_IDX].rank = rank;
384 data_pairs[FRAC_MASK_IDX].len =
386 data_pairs[FRAC_MASK_IDX].rank = rank;
387 data_pairs[COLLECTION_SIZE_IDX].len =
388 ((b == NULL) || (b->collection_size == SIZE_MAX))?
389 -1:(int)b->collection_size;
390 data_pairs[COLLECTION_SIZE_IDX].rank = rank;
391
392 // determine the broadcaster
394 MPI_Allreduce(
395 MPI_IN_PLACE, data_pairs, DATA_PAIR_COUNT, MPI_2INT,
396 MPI_MAXLOC, comm), comm);
397
398 // if at least one process has a valid timestep
399 if (data_pairs[TIMESTEP_IDX].len > 0) {
400
401 // broadcast timestep
402 char * timestep_buffer = xmalloc((size_t)data_pairs[TIMESTEP_IDX].len);
403 if (data_pairs[TIMESTEP_IDX].rank == rank)
404 strcpy(timestep_buffer, b->timestep);
406 MPI_Bcast(
407 timestep_buffer, data_pairs[TIMESTEP_IDX].len, MPI_CHAR,
408 data_pairs[TIMESTEP_IDX].rank, comm), comm);
409
410 // check for consistency
412 (b == NULL) || (b->timestep == NULL) ||
413 !strcmp(b->timestep, timestep_buffer),
414 "ERROR(yac_couple_config_field_merge): "
415 "inconsistent timestep definition (\"%s\" != \"%s\")",
416 b->timestep, timestep_buffer);
417
418 // update timestep
419 free((void*)(a->timestep));
420 a->timestep = timestep_buffer;
421 }
422
423 // if at least one process has a valid metadata
424 if (data_pairs[METADATA_IDX].len > 0) {
425
426 // broadcast metadata
427 char * metadata_buffer = xmalloc((size_t)data_pairs[METADATA_IDX].len);
428 if (data_pairs[METADATA_IDX].rank == rank)
429 strcpy(metadata_buffer, b->metadata);
431 MPI_Bcast(
432 metadata_buffer, data_pairs[METADATA_IDX].len, MPI_CHAR,
433 data_pairs[METADATA_IDX].rank, comm), comm);
434
435 // check for consistency
437 (b == NULL) || (b->metadata == NULL) ||
438 !strcmp(b->metadata, metadata_buffer),
439 "ERROR(yac_couple_config_field_merge): "
440 "inconsistent metadata definition (\"%s\" != \"%s\")",
441 b->metadata, metadata_buffer);
442
443 // update metadata
444 free(a->metadata);
445 a->metadata = metadata_buffer;
446 }
447
448 // if at least one process has a valid fractional mask fallback value
449 if (data_pairs[FRAC_MASK_IDX].len != 0) {
450
451 // broadcast fractional mask fallback value
452 double frac_mask_fallback_value;
453 if (data_pairs[FRAC_MASK_IDX].rank == rank)
454 frac_mask_fallback_value = b->frac_mask_fallback_value;
456 MPI_Bcast(
457 &frac_mask_fallback_value, 1, MPI_DOUBLE,
458 data_pairs[FRAC_MASK_IDX].rank, comm), comm);
459
460 // check for consistency
461 // (use memcmp for comparison, because it can be nan)
463 (b == NULL) || (b->frac_mask_fallback_value == YAC_FRAC_MASK_NO_VALUE) ||
464 !memcmp(
465 &b->frac_mask_fallback_value, &frac_mask_fallback_value,
466 sizeof(double)),
467 "ERROR(yac_couple_config_field_merge): "
468 "inconsistent fractional mask fallback value definition "
469 "(%lf != %lf)",
470 b->frac_mask_fallback_value, frac_mask_fallback_value);
471
472 // update fractional mask fallback value
473 a->frac_mask_fallback_value = frac_mask_fallback_value;
474 }
475
476
477 // if at least one process has a valid collection size
478 if (data_pairs[COLLECTION_SIZE_IDX].len > -1) {
479
480 // check for consistency
482 (b == NULL) || (b->collection_size == SIZE_MAX) ||
483 (b->collection_size == (size_t)(data_pairs[COLLECTION_SIZE_IDX].len)),
484 "ERROR(yac_couple_config_field_merge): "
485 "inconsistent collection size definition (\"%zu\" != \"%d\")",
486 b->collection_size, data_pairs[COLLECTION_SIZE_IDX].len);
487
488 // update collection size
489 a->collection_size = (size_t)(data_pairs[COLLECTION_SIZE_IDX].len);
490 }
491}
492
494 void * a_, void * b_, MPI_Comm comm) {
495
496 UNUSED(comm);
497
498 if (b_ == NULL) return;
499
500 struct yac_couple_config_field_couple * a = a_;
501 struct yac_couple_config_field_couple * b = b_;
502
504 (a->source.lag == b->source.lag),
505 "ERROR(yac_couple_config_field_couple_merge): "
506 "inconsistent source lag (%d != %d)", a->source.lag, b->source.lag)
508 (a->target.lag == b->target.lag),
509 "ERROR(yac_couple_config_field_couple_merge): "
510 "inconsistent target lag (%d != %d)", a->target.lag, b->target.lag)
513 "ERROR(yac_couple_config_field_couple_merge): "
514 "inconsistent mapping side (%d != %d)",
518 "ERROR(yac_couple_config_field_couple_merge): "
519 "inconsistent interpolation stack")
522 "ERROR(yac_couple_config_field_couple_merge): "
523 "inconsistent coupling period operation (%d != %d)",
526 !strcmp(a->coupling_period, b->coupling_period),
527 "ERROR(yac_couple_config_field_couple_merge): "
528 "inconsistent coupling period (%s != %s)",
532 "ERROR(yac_couple_config_field_couple_merge): "
533 "inconsistent defintion of enforce_write_weight_file (%d != %d)",
537 !strcmp(a->weight_file_name, b->weight_file_name),
538 "ERROR(yac_couple_config_field_couple_merge): "
539 "inconsistent weight_file_name (%s != %s)",
542 a->scale_factor == b->scale_factor,
543 "ERROR(yac_couple_config_field_couple_merge): "
544 "inconsistent scale factor (%lf != %lf)",
548 "ERROR(yac_couple_config_field_couple_merge): "
549 "inconsistent scale summand (%lf != %lf)",
553 "ERROR(yac_couple_config_field_couple_merge): "
554 "inconsistent number of source mask names (%zu != %zu)",
557 (a->src_mask_names != NULL) == (b->src_mask_names != NULL),
558 "ERROR(yac_couple_config_field_couple_merge): "
559 "inconsistent availability of source mask names (%s != %s)",
560 (a->src_mask_names != NULL)?"TRUE":"FALSE",
561 (b->src_mask_names != NULL)?"TRUE":"FALSE")
562 for (size_t i = 0; i < a->num_src_mask_names; ++i)
564 !strcmp(a->src_mask_names[i], b->src_mask_names[i]),
565 "ERROR(yac_couple_config_field_couple_merge): "
566 "inconsistent source mask names at index %zu (\"%s\" != \"%s\")",
567 i, a->src_mask_names[i], b->src_mask_names[i])
569 (a->tgt_mask_name != NULL) == (b->tgt_mask_name != NULL),
570 "ERROR(yac_couple_config_field_couple_merge): "
571 "inconsistent availability of target mask name (%s != %s)",
572 (a->tgt_mask_name != NULL)?"TRUE":"FALSE",
573 (b->tgt_mask_name != NULL)?"TRUE":"FALSE")
575 (a->tgt_mask_name == NULL) ||
576 !strcmp(a->tgt_mask_name, b->tgt_mask_name),
577 "ERROR(yac_couple_config_field_couple_merge): "
578 "inconsistent target mask name (\"%s\" != \"%s\")",
580}
581
582static void merge_field_couples(
583 size_t * num_field_couples,
584 struct yac_couple_config_field_couple ** field_couples, MPI_Comm comm);
585
587 void * a_, void * b_, MPI_Comm comm) {
588
589 struct yac_couple_config_couple * a = a_;
590 struct yac_couple_config_couple * b = b_;
591
592 if (b) {
595 b->num_field_couples = 0;
596 b->field_couples = NULL;
597 }
598
599 // distribute and merge field couples
601}
602
603void yac_couple_config_delete(struct yac_couple_config * couple_config) {
604
605 free((void*)couple_config->start_datetime);
606 free((void*)couple_config->end_datetime);
607
608 free((void*)couple_config->config_filename);
609
610 for (size_t i = 0; i < couple_config->num_grids; ++i){
611 yac_couple_config_grid_free(couple_config->grids + i);
612 }
613 free(couple_config->grids);
614
615 for (size_t i = 0;
616 i < couple_config->num_components; ++i)
618 free(couple_config->components);
619
620 for (size_t couple_idx = 0; couple_idx < couple_config->num_couples;
621 ++couple_idx)
622 yac_couple_config_couple_free(couple_config->couples + couple_idx);
623 free(couple_config->couples);
624 free(couple_config);
625}
626
628 struct yac_couple_config * couple_config, char const * name) {
629
631 name != NULL,
632 "ERROR(yac_couple_config_add_grid_): "
633 "invalid name (NULL is not allowed)")
634
635 for (size_t i = 0; i < couple_config->num_grids; ++i)
636 if (!strcmp(couple_config->grids[i].name, name)) return i;
637
638 size_t grid_idx = couple_config->num_grids;
639 couple_config->num_grids++;
640 couple_config->grids =
641 xrealloc(
642 couple_config->grids,
643 couple_config->num_grids * sizeof(*(couple_config->grids)));
644 couple_config->grids[grid_idx].name = string_dup(name);
645 couple_config->grids[grid_idx].metadata = NULL;
646 return grid_idx;
647}
648
650 struct yac_couple_config * couple_config, char const * name) {
651
652 yac_couple_config_add_grid_(couple_config, name);
653}
654
656 struct yac_couple_config * couple_config, char const * name) {
657
659 name != NULL,
660 "ERROR(yac_couple_config_add_component_): "
661 "invalid name (NULL is not allowed)")
662
663 for (size_t i = 0; i < couple_config->num_components; ++i)
664 if (!strcmp(couple_config->components[i].name, name)) return i;
665
666 size_t component_idx = couple_config->num_components;
667 couple_config->num_components++;
668 couple_config->components =
669 xrealloc(
670 couple_config->components,
671 couple_config->num_components * sizeof(*(couple_config->components)));
672 couple_config->components[component_idx].fields = NULL;
673 couple_config->components[component_idx].num_fields = 0;
674 couple_config->components[component_idx].name = strdup(name);
675 couple_config->components[component_idx].metadata = NULL;
676 return component_idx;
677}
678
680 struct yac_couple_config * couple_config, char const * name) {
681
682 yac_couple_config_add_component_(couple_config, name);
683}
684
686 struct yac_couple_config * couple_config,
687 char const * comp_name, const char* metadata) {
688 size_t comp_idx = yac_couple_config_get_component_idx(couple_config, comp_name);
689 free(couple_config->components[comp_idx].metadata);
690 couple_config->components[comp_idx].metadata
691 = metadata==NULL?NULL:strdup(metadata);
692}
693
695 struct yac_couple_config * couple_config,
696 char const * grid_name, const char* metadata) {
697 if(!yac_couple_config_contains_grid_name(couple_config, grid_name))
698 yac_couple_config_add_grid(couple_config, grid_name);
699 size_t grid_idx = yac_couple_config_get_grid_idx(couple_config, grid_name);
700 free(couple_config->grids[grid_idx].metadata);
701 couple_config->grids[grid_idx].metadata
702 = metadata==NULL?NULL:strdup(metadata);
703}
704
706 struct yac_couple_config * couple_config,
707 const char* comp_name, const char * grid_name, const char* field_name,
708 const char* metadata) {
709 size_t comp_idx = yac_couple_config_get_component_idx(couple_config, comp_name);
710 size_t grid_idx = yac_couple_config_get_grid_idx(couple_config, grid_name);
711 size_t field_idx = yac_couple_config_get_field_idx(couple_config,
712 comp_idx, grid_idx, field_name);
713 free(couple_config->components[comp_idx].fields[field_idx].metadata);
714 couple_config->components[comp_idx].fields[field_idx].metadata
715 = metadata==NULL?NULL:strdup(metadata);
716}
717
719 struct yac_couple_config * couple_config,
720 const char * comp_name) {
721 size_t comp_idx = yac_couple_config_get_component_idx(couple_config, comp_name);
722 return couple_config->components[comp_idx].metadata;
723}
724
726 struct yac_couple_config * couple_config,
727 const char * grid_name) {
728 size_t grid_idx = yac_couple_config_get_grid_idx(couple_config, grid_name);
729 return couple_config->grids[grid_idx].metadata;
730}
731
733 struct yac_couple_config * couple_config,
734 const char* comp_name, const char * grid_name, const char* field_name) {
735 size_t comp_idx = yac_couple_config_get_component_idx(couple_config, comp_name);
736 size_t grid_idx = yac_couple_config_get_grid_idx(couple_config, grid_name);
737 size_t field_idx = yac_couple_config_get_field_idx(couple_config,
738 comp_idx, grid_idx, field_name);
739 return couple_config->components[comp_idx].fields[field_idx].metadata;
740}
741
743 struct yac_couple_config * couple_config, size_t component_idx,
744 char const * routine_name, int line) {
745
747 component_idx < couple_config->num_components,
748 "ERROR(%s:%d:%s): invalid component_idx", __FILE__, line, routine_name)
749}
750
751static void check_grid_idx(
752 struct yac_couple_config * couple_config, size_t grid_idx,
753 char const * routine_name, int line) {
754
756 grid_idx < couple_config->num_grids,
757 "ERROR(%s:%d:%s): invalid grid_idx", __FILE__, line, routine_name)
758}
759
761 struct yac_couple_config * couple_config,
762 size_t comp_idx, size_t grid_idx, char const * name,
763 char const * timestep, size_t collection_size) {
764
766 couple_config, comp_idx,
767 "yac_couple_config_component_add_field_", __LINE__);
769 couple_config, grid_idx,
770 "yac_couple_config_component_add_field_", __LINE__);
771
772 // check whether the field already exists
773 struct yac_couple_config_component * component =
774 couple_config->components + comp_idx;
775 for (size_t i = 0; i < component->num_fields; i++) {
776 if((strcmp(component->fields[i].name, name) == 0) &&
777 (component->fields[i].grid_idx == grid_idx)) {
778 // if no timestep is defined for the field
779 if (!component->fields[i].timestep && timestep)
780 component->fields[i].timestep = strdup(timestep);
781 // if no collection size is defined for the field
782 if (component->fields[i].collection_size == SIZE_MAX)
783 component->fields[i].collection_size = collection_size;
785 !timestep ||
786 !strcmp(timestep, component->fields[i].timestep),
787 "ERROR(yac_couple_config_component_add_field): "
788 "inconsistent timestep definition (\"%s\" != \"%s\")",
789 timestep, component->fields[i].timestep);
791 collection_size == SIZE_MAX ||
792 collection_size == component->fields[i].collection_size,
793 "ERROR(yac_couple_config_component_add_field): "
794 "inconsistent collection_size definition (%zu != %zu)",
795 collection_size, component->fields[i].collection_size);
796 return i;
797 }
798 }
799
800 size_t field_idx = component->num_fields;
801 component->num_fields++;
802
803 component->fields =
804 xrealloc(
805 component->fields,
806 component->num_fields *
807 sizeof(*(component->fields)));
808 struct yac_couple_config_field * field =
809 component->fields + field_idx;
810
811 field->name = strdup(name);
812 field->grid_idx = grid_idx;
813 field->timestep = timestep?strdup(timestep):NULL;
814 field->metadata = NULL;
817 return field_idx;
818}
819
821 struct yac_couple_config * couple_config, const char* component_name,
822 const char* grid_name, const char* name, char const * timestep,
823 size_t collection_size) {
824
826 couple_config,
827 yac_couple_config_get_component_idx(couple_config, component_name),
828 yac_couple_config_get_grid_idx(couple_config, grid_name),
830}
831
833 struct yac_couple_config * couple_config) {
834
835 return couple_config->num_couples;
836}
837
839 struct yac_couple_config * couple_config, size_t couple_idx,
840 char const * routine_name, int line) {
841
843 couple_idx < couple_config->num_couples,
844 "ERROR(%s:%d:%s): invalid couple_idx", __FILE__, line, routine_name);
845}
846
848 struct yac_couple_config * couple_config, size_t couple_idx) {
849
851 couple_config, couple_idx, "yac_couple_config_get_num_couple_fields",
852 __LINE__);
853
854 return couple_config->couples[couple_idx].num_field_couples;
855}
856
858 struct yac_couple_config * couple_config, size_t couple_idx,
859 char const * couple_component_names[2]) {
860
862 couple_config, couple_idx, "yac_couple_config_get_couple_component_names",
863 __LINE__);
864
865 for (int i = 0; i < 2; ++i)
866 couple_component_names[i] =
867 couple_config->components[
868 couple_config->couples[couple_idx].component_indices[i]].name;
869}
870
872 struct yac_couple_config * couple_config, char const * component_name) {
873
875 component_name,
876 "ERROR(yac_couple_config_component_name_is_valid): component name is NULL")
878 strlen(component_name) <= YAC_MAX_CHARLEN,
879 "ERROR(yac_couple_config_component_name_is_valid): "
880 "component name is too long (maximum is YAC_MAX_CHARLEN)")
881
882 for (size_t component_idx = 0; component_idx < couple_config->num_components;
883 ++component_idx)
884 if (!strcmp(component_name, couple_config->components[component_idx].name))
885 return 1;
886 return 0;
887}
888
890 struct yac_couple_config * couple_config) {
891
892 return couple_config->num_components;
893}
894
896 struct yac_couple_config * couple_config) {
897
898 return couple_config->num_grids;
899}
900
902 struct yac_couple_config * couple_config, size_t component_idx) {
903 check_component_idx(couple_config, component_idx,
904 "yac_couple_config_get_num_fields", __LINE__);
905 return couple_config->components[component_idx].num_fields;
906}
907
909 struct yac_couple_config * couple_config, char const * component_name) {
910
911 size_t component_idx = SIZE_MAX;
912 for (size_t i = 0; (i < couple_config->num_components) &&
913 (component_idx == SIZE_MAX); ++i)
914 if (!strcmp(couple_config->components[i].name, component_name))
915 component_idx = i;
916
918 component_idx != SIZE_MAX,
919 "ERROR(yac_couple_config_get_component_idx): "
920 "Component \"%s\" not found in coupling configuration",
921 component_name);
922
923 return component_idx;
924}
925
927 struct yac_couple_config * couple_config, char const * grid_name) {
928
929 size_t grid_idx = SIZE_MAX;
930 for (size_t i = 0;
931 (i < couple_config->num_grids) && (grid_idx == SIZE_MAX); ++i)
932 if (!strcmp(couple_config->grids[i].name, grid_name))
933 grid_idx = i;
934
936 grid_idx != SIZE_MAX,
937 "ERROR(yac_couple_config_get_grid_idx): "
938 "grid name \"%s\" not in list of grids", grid_name)
939
940 return grid_idx;
941}
942
944 struct yac_couple_config * couple_config, size_t component_idx,
945 size_t grid_idx, char const * field_name) {
947 couple_config, component_idx,
948 "yac_couple_config_get_component_name", __LINE__);
949
950 size_t field_idx = SIZE_MAX;
951 struct yac_couple_config_component * component =
952 couple_config->components + component_idx;
953 size_t nbr_fields = component->num_fields;
954 for(size_t i=0;
955 (i<nbr_fields) && (field_idx == SIZE_MAX); ++i)
956 if((component->fields[i].grid_idx == grid_idx) &&
957 !strcmp(
958 field_name,
959 component->fields[i].name))
960 field_idx = i;
961
963 field_idx != SIZE_MAX,
964 "ERROR(yac_couple_config_get_field_idx): "
965 "field not found "
966 "(component_idx %zu grid_idx %zu field_name \"%s\"",
967 component_idx, grid_idx, field_name);
968
969 return field_idx;
970}
971
973 struct yac_couple_config * couple_config, size_t component_idx) {
974
976 couple_config, component_idx,
977 "yac_couple_config_get_component_name", __LINE__);
978
979 return couple_config->components[component_idx].name;
980}
981
982static void check_field_idx(
983 struct yac_couple_config * couple_config, size_t component_idx,
984 size_t field_idx, char const * routine_name, int line) {
985
987 couple_config, component_idx, routine_name, __LINE__);
988
990 field_idx <
991 couple_config->components[component_idx].num_fields,
992 "ERROR(%s:%d:%s): invalid field_idx", __FILE__,
993 line, routine_name)
994}
995
997 struct yac_couple_config * couple_config, size_t component_idx,
998 size_t field_idx) {
999
1001 couple_config, component_idx, field_idx,
1002 "yac_couple_config_get_field_grid_name", __LINE__);
1003
1004 return
1005 couple_config->grids[
1006 couple_config->
1007 components[component_idx].
1008 fields[field_idx].grid_idx].name;
1009}
1010
1012 struct yac_couple_config * couple_config, size_t component_idx,
1013 size_t field_idx) {
1014
1016 couple_config, component_idx, field_idx,
1017 "yac_couple_config_get_field_name", __LINE__);
1018
1019 return
1020 couple_config->
1021 components[component_idx].
1022 fields[field_idx].name;
1023}
1024
1026 struct yac_couple_config * couple_config,
1027 char const * component_name, char const * grid_name,
1028 char const * field_name) {
1029
1030 size_t component_idx =
1031 yac_couple_config_get_component_idx(couple_config, component_name);
1032 size_t grid_idx =
1033 yac_couple_config_get_grid_idx(couple_config, grid_name);
1034 size_t field_idx =
1036 couple_config, component_idx, grid_idx, field_name);
1037
1038 struct yac_couple_config_field * field =
1039 couple_config->components[component_idx].fields +
1040 field_idx;
1041
1043 field->timestep,
1044 "ERROR(yac_couple_config_get_field_timestep): "
1045 "no valid timestep defined (component: \"%s\" field \"%s\")",
1046 couple_config->components[component_idx].name, field->name);
1047
1048 return field->timestep;
1049}
1050
1052 struct yac_couple_config * couple_config,
1053 char const * component_name, char const * grid_name,
1054 char const * field_name) {
1055
1056 size_t component_idx =
1057 yac_couple_config_get_component_idx(couple_config, component_name);
1058 size_t grid_idx =
1059 yac_couple_config_get_grid_idx(couple_config, grid_name);
1060 size_t field_idx =
1062 couple_config, component_idx, grid_idx, field_name);
1063
1064 size_t nbr_couples = couple_config->num_couples;
1065 for(size_t couple_idx = 0; couple_idx<nbr_couples; ++couple_idx){
1066 struct yac_couple_config_couple * couple = couple_config->couples + couple_idx;
1067 if(couple->component_indices[0] != component_idx &&
1068 couple->component_indices[1] != component_idx)
1069 continue;
1070 size_t nbr_trans_couples = couple->num_field_couples;
1071 for(size_t trans_couple_idx = 0; trans_couple_idx < nbr_trans_couples;
1072 ++trans_couple_idx){
1073 struct yac_couple_config_field_couple * transcouple =
1074 couple->field_couples + trans_couple_idx;
1075 if(transcouple->source.component_idx == component_idx &&
1076 transcouple->source.field_idx == field_idx)
1078 if(transcouple->target.component_idx == component_idx &&
1079 transcouple->target.field_idx == field_idx)
1081 }
1082 }
1084}
1085
1087 struct yac_couple_config * couple_config,
1088 size_t component_idx, size_t field_idx) {
1089
1091 couple_config, component_idx, field_idx,
1092 "yac_couple_config_field_is_valid", __LINE__);
1093
1094 struct yac_couple_config_field * field =
1095 couple_config->components[component_idx].fields +
1096 field_idx;
1097
1098 return (field->collection_size != SIZE_MAX) &&
1099 (field->timestep != NULL);
1100}
1101
1103 struct yac_couple_config * couple_config, size_t couple_idx,
1104 size_t field_couple_idx, char const * routine_name, int line) {
1105
1106 check_couple_idx(couple_config, couple_idx, routine_name, __LINE__);
1107
1109 field_couple_idx <
1110 couple_config->couples[couple_idx].num_field_couples,
1111 "ERROR(%s:%d:%s): invalid field_couple_idx",
1112 __FILE__, line, routine_name)
1113}
1114
1116 struct yac_couple_config * couple_config,
1117 size_t couple_idx, size_t field_couple_idx) {
1118
1120 couple_config, couple_idx, field_couple_idx,
1121 "yac_couple_config_get_interp_stack", __LINE__);
1122
1123 return
1124 couple_config->
1125 couples[couple_idx].
1126 field_couples[field_couple_idx].interp_stack;
1127}
1128
1130 struct yac_couple_config * couple_config,
1131 size_t couple_idx, size_t field_couple_idx,
1132 char const ** src_grid_name, char const ** tgt_grid_name) {
1133
1135 couple_config, couple_idx, field_couple_idx,
1136 "yac_couple_config_get_field_grid_names", __LINE__);
1137
1138 size_t src_component_idx =
1139 couple_config->couples[couple_idx].
1140 field_couples[field_couple_idx].
1141 source.component_idx;
1142 size_t src_field_idx =
1143 couple_config->couples[couple_idx].
1144 field_couples[field_couple_idx].
1145 source.field_idx;
1146
1147 size_t tgt_component_idx =
1148 couple_config->couples[couple_idx].
1149 field_couples[field_couple_idx].
1150 target.component_idx;
1151 size_t tgt_field_idx =
1152 couple_config->couples[couple_idx].
1153 field_couples[field_couple_idx].
1154 target.field_idx;
1155
1156 *src_grid_name =
1157 couple_config->grids[
1158 couple_config->components[src_component_idx].
1159 fields[src_field_idx].grid_idx].name;
1160 *tgt_grid_name =
1161 couple_config->grids[
1162 couple_config->components[tgt_component_idx].
1163 fields[tgt_field_idx].grid_idx].name;
1164}
1165
1167 struct yac_couple_config * couple_config,
1168 char const * comp_name, char const * grid_name, char const * field_name,
1169 double frac_mask_fallback_value) {
1170
1172 (frac_mask_fallback_value != YAC_FRAC_MASK_UNDEF) &&
1173 (frac_mask_fallback_value != YAC_FRAC_MASK_NO_VALUE),
1174 "ERROR(yac_couple_config_field_enable_frac_mask): "
1175 "\"%lf\" is not a valid fractional mask fallback value "
1176 "(component: \"%s\" grid: \"%s\" field \"%s\")",
1177 frac_mask_fallback_value, comp_name, grid_name, field_name);
1178
1179 size_t comp_idx = yac_couple_config_get_component_idx(couple_config, comp_name);
1180 size_t grid_idx = yac_couple_config_get_grid_idx(couple_config, grid_name);
1181 size_t field_idx = yac_couple_config_get_field_idx(couple_config,
1182 comp_idx, grid_idx, field_name);
1183
1184 double old_frac_mask_fallback_value =
1185 couple_config->components[comp_idx].fields[field_idx].
1186 frac_mask_fallback_value;
1187
1189 (old_frac_mask_fallback_value == YAC_FRAC_MASK_UNDEF) ||
1190 (old_frac_mask_fallback_value == YAC_FRAC_MASK_NO_VALUE) ||
1191 (old_frac_mask_fallback_value == frac_mask_fallback_value),
1192 "ERROR(yac_couple_config_field_enable_frac_mask): "
1193 "fractional mask fallback value was already set:\n"
1194 "\told value \"%lf\" new value \"%lf\" "
1195 "(component: \"%s\" grid: \"%s\" field \"%s\")",
1196 old_frac_mask_fallback_value, frac_mask_fallback_value,
1197 comp_name, grid_name, field_name);
1198
1199 couple_config->components[comp_idx].fields[field_idx].
1200 frac_mask_fallback_value = frac_mask_fallback_value;
1201}
1202
1204 struct yac_couple_config * couple_config,
1205 char const * component_name, char const * grid_name,
1206 char const * field_name) {
1207
1208 size_t component_idx =
1209 yac_couple_config_get_component_idx(couple_config, component_name);
1210 size_t grid_idx =
1211 yac_couple_config_get_grid_idx(couple_config, grid_name);
1212 size_t field_idx =
1214 couple_config, component_idx, grid_idx, field_name);
1215
1216 struct yac_couple_config_field * field =
1217 couple_config->components[component_idx].fields +
1218 field_idx;
1219
1222 "ERROR(yac_couple_config_get_frac_mask_fallback_value): "
1223 "no valid fractional mask fallback value defined "
1224 "(component: \"%s\" grid: \"%s\" field \"%s\")",
1225 couple_config->components[component_idx].name, grid_name, field->name);
1226
1227 return field->frac_mask_fallback_value;
1228}
1229
1231 struct yac_couple_config * couple_config,
1232 char const * component_name, char const * grid_name,
1233 char const * field_name) {
1234
1235 size_t component_idx =
1236 yac_couple_config_get_component_idx(couple_config, component_name);
1237 size_t grid_idx =
1238 yac_couple_config_get_grid_idx(couple_config, grid_name);
1239 size_t field_idx =
1241 couple_config, component_idx, grid_idx, field_name);
1242
1243 struct yac_couple_config_field * field =
1244 couple_config->components[component_idx].fields +
1245 field_idx;
1246
1248 field->collection_size != SIZE_MAX,
1249 "ERROR(yac_couple_config_get_collection_size): "
1250 "no valid collection size defined "
1251 "(component: \"%s\" grid: \"%s\" field \"%s\")",
1252 couple_config->components[component_idx].name, grid_name, field->name);
1253
1254 return field->collection_size;
1255}
1256
1258 struct yac_couple_config * couple_config,
1259 size_t couple_idx, size_t field_couple_idx,
1260 char const ** src_component_name, char const ** tgt_component_name) {
1261
1263 couple_config, couple_idx, field_couple_idx,
1264 "yac_couple_config_get_field_couple_component_names", __LINE__);
1265
1266 *src_component_name =
1267 couple_config->components[
1268 couple_config->couples[couple_idx].
1269 field_couples[field_couple_idx].source.component_idx].name;
1270 *tgt_component_name =
1271 couple_config->components[
1272 couple_config->couples[couple_idx].
1273 field_couples[field_couple_idx].target.component_idx].name;
1274}
1275
1277 struct yac_couple_config * couple_config,
1278 size_t couple_idx, size_t field_couple_idx,
1279 char const ** src_field_name, const char ** tgt_field_name) {
1280
1282 couple_config, couple_idx, field_couple_idx,
1283 "yac_couple_config_get_field_names", __LINE__);
1284
1285 size_t src_component_idx =
1286 couple_config->couples[couple_idx].
1287 field_couples[field_couple_idx].source.component_idx;
1288 size_t tgt_component_idx =
1289 couple_config->couples[couple_idx].
1290 field_couples[field_couple_idx].target.component_idx;
1291 size_t src_field_idx =
1292 couple_config->couples[couple_idx].
1293 field_couples[field_couple_idx].source.field_idx;
1294 size_t tgt_field_idx =
1295 couple_config->couples[couple_idx].
1296 field_couples[field_couple_idx].target.field_idx;
1297
1298 *src_field_name =
1299 couple_config->components[src_component_idx].
1300 fields[src_field_idx].name;
1301 *tgt_field_name =
1302 couple_config->components[tgt_component_idx].
1303 fields[tgt_field_idx].name;
1304}
1305
1307 struct yac_couple_config * couple_config,
1308 size_t couple_idx, size_t field_couple_idx) {
1309
1311 couple_config, couple_idx, field_couple_idx,
1312 "yac_couple_config_mapping_on_source", __LINE__);
1313
1314 return
1315 couple_config->
1316 couples[couple_idx].
1317 field_couples[field_couple_idx].
1318 mapping_on_source;
1319}
1320
1322 struct yac_couple_config * couple_config,
1323 size_t couple_idx, size_t field_couple_idx) {
1324
1326 couple_config, couple_idx, field_couple_idx,
1327 "yac_couple_config_get_source_lag", __LINE__);
1328
1329 return
1330 couple_config->
1331 couples[couple_idx].
1332 field_couples[field_couple_idx].
1333 source.lag;
1334}
1335
1337 struct yac_couple_config * couple_config,
1338 size_t couple_idx, size_t field_couple_idx) {
1339
1341 couple_config, couple_idx, field_couple_idx,
1342 "yac_couple_config_get_target_lag", __LINE__);
1343
1344 return
1345 couple_config->
1346 couples[couple_idx].
1347 field_couples[field_couple_idx].
1348 target.lag;
1349}
1350
1352 struct yac_couple_config * couple_config,
1353 size_t couple_idx, size_t field_couple_idx) {
1354
1356 couple_config, couple_idx, field_couple_idx,
1357 "yac_couple_config_get_coupling_period", __LINE__);
1358
1359 return
1360 couple_config->
1361 couples[couple_idx].
1362 field_couples[field_couple_idx].
1363 coupling_period;
1364}
1365
1367 struct yac_couple_config * couple_config,
1368 size_t couple_idx, size_t field_couple_idx) {
1369
1371 couple_config, couple_idx, field_couple_idx,
1372 "yac_couple_config_get_source_timestep", __LINE__);
1373
1374 struct yac_couple_config_field_couple * tcouple =
1375 couple_config->couples[couple_idx].field_couples + field_couple_idx;
1376 char const * timestep =
1377 couple_config->components[tcouple->source.component_idx].
1378 fields[tcouple->source.field_idx].timestep;
1379
1381 timestep,
1382 "ERROR(yac_couple_config_get_source_timestep): "
1383 "no valid timestep defined (component: \"%s\" field \"%s\")",
1384 couple_config->components[tcouple->source.component_idx].name,
1385 couple_config->components[tcouple->source.component_idx].
1386 fields[tcouple->source.field_idx].name);
1387
1388 return timestep;
1389}
1390
1392 struct yac_couple_config * couple_config,
1393 size_t couple_idx, size_t field_couple_idx) {
1394
1396 couple_config, couple_idx, field_couple_idx,
1397 "yac_couple_config_get_target_timestep", __LINE__);
1398
1399 struct yac_couple_config_field_couple * tcouple =
1400 couple_config->couples[couple_idx].field_couples + field_couple_idx;
1401 char const * timestep =
1402 couple_config->components[tcouple->target.component_idx].
1403 fields[tcouple->target.field_idx].timestep;
1404
1406 timestep,
1407 "ERROR(yac_couple_config_get_target_timestep): "
1408 "no valid timestep defined (component: \"%s\" field \"%s\")",
1409 couple_config->components[tcouple->target.component_idx].name,
1410 couple_config->components[tcouple->target.component_idx].
1411 fields[tcouple->target.field_idx].name);
1412
1413 return timestep;
1414}
1415
1417 struct yac_couple_config * couple_config,
1418 size_t couple_idx, size_t field_couple_idx) {
1419
1421 couple_config, couple_idx, field_couple_idx,
1422 "yac_couple_config_get_coupling_period_operation", __LINE__);
1423
1424 return
1425 couple_config->
1426 couples[couple_idx].
1427 field_couples[field_couple_idx].
1429}
1430
1432 struct yac_couple_config * couple_config,
1433 char const * start, char const * end) {
1434
1435 if ((start != NULL) && (strlen(start) > 0)) {
1436
1437 // in case start points to couple_config->start_datetime,
1438 // we have to first strdup start, before freeing it
1439 char const * old_start = couple_config->start_datetime;
1440 couple_config->start_datetime = strdup(start);
1441 free((void*)old_start);
1442 }
1443
1444 if ((end != NULL) && (strlen(end) > 0)) {
1445
1446 // in case end points to couple_config->end_datetime,
1447 // we have to first strdup end, before freeing it
1448 char const * old_end = couple_config->end_datetime;
1449 couple_config->end_datetime = strdup(end);
1450 free((void*)old_end);
1451 }
1452}
1453
1455 struct yac_couple_config * couple_config) {
1456
1457 YAC_ASSERT(couple_config->start_datetime,
1458 "ERROR(yac_couple_config_get_start_datetime): "
1459 "start_datetime not yet defined");
1460
1461 return couple_config->start_datetime;
1462}
1463
1465 struct yac_couple_config * couple_config) {
1466
1467 YAC_ASSERT(couple_config->end_datetime,
1468 "ERROR(yac_couple_config_get_start_datetime): "
1469 "start_datetime not yet defined");
1470
1471 return couple_config->end_datetime;
1472}
1473
1475 struct yac_couple_config * couple_config, size_t grid_idx) {
1476
1478 grid_idx < couple_config->num_grids,
1479 "ERROR(yac_couple_config_get_grid_name): "
1480 "Invalid grid idx %zu", grid_idx);
1481
1482 return couple_config->grids[grid_idx].name;
1483}
1484
1486 struct yac_couple_config * couple_config,
1487 size_t couple_idx, size_t field_couple_idx) {
1488
1490 couple_config, couple_idx, field_couple_idx,
1491 "yac_couple_config_enforce_write_weight_file", __LINE__);
1492
1493 return
1494 couple_config->
1495 couples[couple_idx].
1496 field_couples[field_couple_idx].
1498}
1499
1501 struct yac_couple_config * couple_config,
1502 size_t couple_idx, size_t field_couple_idx) {
1503
1505 couple_config, couple_idx, field_couple_idx,
1506 "yac_couple_config_get_weight_file_name", __LINE__);
1507
1508 static char dummy[] = "\0";
1509 char const * weight_file_name =
1510 couple_config->
1511 couples[couple_idx].
1512 field_couples[field_couple_idx].
1514
1515 return (weight_file_name != NULL)?weight_file_name:dummy;
1516}
1517
1519 struct yac_couple_config * couple_config,
1520 size_t couple_idx, size_t field_couple_idx) {
1521
1523 couple_config, couple_idx, field_couple_idx,
1524 "yac_couple_config_get_scale_factor", __LINE__);
1525
1526 return
1527 couple_config->
1528 couples[couple_idx].
1529 field_couples[field_couple_idx].
1531}
1532
1534 struct yac_couple_config * couple_config,
1535 size_t couple_idx, size_t field_couple_idx) {
1536
1538 couple_config, couple_idx, field_couple_idx,
1539 "yac_couple_config_get_scale_summand", __LINE__);
1540
1541 return
1542 couple_config->
1543 couples[couple_idx].
1544 field_couples[field_couple_idx].
1546}
1547
1549 struct yac_couple_config * couple_config,
1550 size_t couple_idx, size_t field_couple_idx,
1551 char const * const ** mask_names, size_t * num_mask_names) {
1552
1554 couple_config, couple_idx, field_couple_idx,
1555 "yac_couple_config_get_src_mask_names", __LINE__);
1556
1557 *mask_names =
1558 (char const * const *)(
1559 couple_config->
1560 couples[couple_idx].
1561 field_couples[field_couple_idx].
1563 *num_mask_names =
1564 couple_config->
1565 couples[couple_idx].
1566 field_couples[field_couple_idx].
1568}
1569
1571 struct yac_couple_config * couple_config,
1572 size_t couple_idx, size_t field_couple_idx) {
1573
1575 couple_config, couple_idx, field_couple_idx,
1576 "yac_couple_config_get_tgt_mask_name", __LINE__);
1577
1578 return
1579 couple_config->
1580 couples[couple_idx].
1581 field_couples[field_couple_idx].
1583}
1584
1586 struct yac_couple_config * couple_config, char const * grid_name) {
1587
1588 for (size_t grid_idx = 0; grid_idx < couple_config->num_grids;
1589 ++grid_idx)
1590 if (!strcmp(couple_config->grids[grid_idx].name, grid_name))
1591 return 1;
1592 return 0;
1593}
1594
1596 char const * string, MPI_Comm comm) {
1597
1598 int strlen_pack_size, string_pack_size;
1599 yac_mpi_call(MPI_Pack_size(1, MPI_INT, comm, &strlen_pack_size), comm);
1600
1601 if (string != NULL) {
1603 MPI_Pack_size(
1604 (int)(strlen(string)), MPI_CHAR, comm, &string_pack_size), comm);
1605 } else {
1606 string_pack_size = 0;
1607 }
1608
1609 return (size_t)strlen_pack_size + (size_t)string_pack_size;
1610}
1611
1613 void * grid_, MPI_Comm comm) {
1614
1615 struct yac_couple_config_grid * grid = grid_;
1616
1617 return yac_couple_config_get_string_pack_size(grid->name, comm);
1618}
1619
1621 size_t num_grids, void * grids_, MPI_Comm comm) {
1622
1623 struct yac_couple_config_grid * grids = grids_;
1624
1625 int num_grids_pack_size;
1626 yac_mpi_call(MPI_Pack_size(1, MPI_INT, comm, &num_grids_pack_size), comm);
1627
1628 size_t grids_pack_size = 0;
1629 for (size_t i = 0; i < num_grids; ++i)
1630 grids_pack_size +=
1632
1633 return (size_t)num_grids_pack_size + grids_pack_size;
1634}
1635
1637 struct yac_couple_config_field * field,
1638 MPI_Comm comm) {
1639
1640 YAC_ASSERT(
1641 field->grid_idx <= INT_MAX,
1642 "ERROR(yac_couple_config_get_field_pack_size):"
1643 "grid_idx is too big")
1644
1645 int int_pack_size;
1646 yac_mpi_call(MPI_Pack_size(1, MPI_INT, comm, &int_pack_size), comm);
1647 int double_pack_size;
1648 yac_mpi_call(MPI_Pack_size(1, MPI_DOUBLE, comm, &double_pack_size), comm);
1649 size_t name_pack_size = yac_couple_config_get_string_pack_size(
1650 field->name, comm);
1651 size_t timestep_pack_size = yac_couple_config_get_string_pack_size(
1652 field->timestep, comm);
1653 size_t metadata_pack_size = yac_couple_config_get_string_pack_size(
1654 field->metadata, comm);
1655
1656 return int_pack_size + // grid_idx
1657 int_pack_size + // collection_size
1658 double_pack_size + // frac_mask_fallback_value
1659 name_pack_size +
1660 timestep_pack_size +
1661 metadata_pack_size;
1662}
1663
1665 size_t num_fields, void * fields_, MPI_Comm comm) {
1666
1667 struct yac_couple_config_field * fields = fields_;
1668
1669 int num_fields_pack_size;
1671 MPI_Pack_size(1, MPI_INT, comm, &num_fields_pack_size), comm);
1672
1673 size_t fields_pack_size = 0;
1674 for (size_t i = 0; i < num_fields; ++i)
1675 fields_pack_size +=
1677 fields + i, comm);
1678
1679 return (size_t)num_fields_pack_size +
1680 fields_pack_size;
1681}
1682
1684 struct yac_couple_config_component * component, MPI_Comm comm) {
1685
1686 return
1688}
1689
1691 size_t num_components, void * components_, MPI_Comm comm) {
1692
1693 struct yac_couple_config_component * components = components_;
1694
1695 int num_components_pack_size;
1697 MPI_Pack_size(1, MPI_INT, comm, &num_components_pack_size), comm);
1698
1699 size_t components_pack_size = 0;
1700 for (size_t i = 0; i < num_components; ++i)
1701 components_pack_size +=
1703
1704 return (size_t)num_components_pack_size + components_pack_size;
1705}
1706
1708 struct yac_couple_config_field_couple * field_couple, MPI_Comm comm) {
1709
1710 int ints_pack_size;
1711 yac_mpi_call(MPI_Pack_size(10, MPI_INT, comm, &ints_pack_size), comm);
1712 int doubles_pack_size;
1713 yac_mpi_call(MPI_Pack_size(2, MPI_DOUBLE, comm, &doubles_pack_size), comm);
1714 int src_mask_names_pack_size = 0;
1715 if (field_couple->num_src_mask_names > 0)
1716 for (size_t i = 0; i < field_couple->num_src_mask_names; ++i)
1717 src_mask_names_pack_size +=
1719 field_couple->src_mask_names[i], comm);
1720
1721 return
1722 (size_t)ints_pack_size + // source.component_idx
1723 // source.field_idx
1724 // source.lag
1725 // target.component_idx
1726 // target.field_idx
1727 // target.lag
1728 // mapping_on_source
1729 // coupling_period_operation
1730 // enforce_write_weight_file
1731 // num_src_mask_names
1733 field_couple->interp_stack, comm) +
1735 field_couple->coupling_period, comm) +
1737 field_couple->weight_file_name, comm) +
1738 doubles_pack_size + // scale_factor
1739 // scale_summand
1740 src_mask_names_pack_size +
1742 field_couple->tgt_mask_name, comm);
1743}
1744
1746 size_t num_field_couples, void * field_couples_, MPI_Comm comm) {
1747
1748 struct yac_couple_config_field_couple * field_couples = field_couples_;
1749
1750 int num_field_couples_pack_size;
1752 MPI_Pack_size(1, MPI_INT, comm, &num_field_couples_pack_size), comm);
1753
1754 size_t field_couples_pack_size = 0;
1755 for (size_t i = 0; i < num_field_couples; ++i)
1756 field_couples_pack_size +=
1757 yac_couple_config_get_field_couple_pack_size(field_couples + i, comm);
1758
1759 return (size_t)num_field_couples_pack_size +
1760 field_couples_pack_size;
1761}
1762
1764 struct yac_couple_config_couple * couple, MPI_Comm comm) {
1765
1766 UNUSED(couple);
1767
1768 int component_indices_pack_size;
1770 MPI_Pack_size(2, MPI_INT, comm, &component_indices_pack_size), comm);
1771
1772 return (size_t)component_indices_pack_size;
1773}
1774
1776 size_t num_couples, void * couples_, MPI_Comm comm) {
1777
1778 struct yac_couple_config_couple * couples = couples_;
1779
1780 int num_couples_pack_size;
1782 MPI_Pack_size(1, MPI_INT, comm, &num_couples_pack_size), comm);
1783
1784 size_t couples_pack_size = 0;
1785 for (size_t i = 0; i < num_couples; ++i)
1786 couples_pack_size +=
1788
1789 return (size_t)num_couples_pack_size + couples_pack_size;
1790}
1791
1793 char const * string, void * buffer, int buffer_size, int * position,
1794 MPI_Comm comm) {
1795
1796 size_t len = (string == NULL)?0:strlen(string);
1797
1798 YAC_ASSERT(
1799 len <= INT_MAX, "ERROR(yac_couple_config_pack_string): string too long")
1800
1801 int len_int = (int)len;
1802
1804 MPI_Pack(
1805 &len_int, 1, MPI_INT, buffer, buffer_size, position, comm), comm);
1806
1807 if (len > 0)
1809 MPI_Pack(
1810 string, len_int, MPI_CHAR, buffer, buffer_size, position, comm),
1811 comm);
1812}
1813
1815 struct yac_couple_config_grid * grid,
1816 void * buffer, int buffer_size, int * position, MPI_Comm comm) {
1817
1819 grid->name, buffer, buffer_size, position, comm);
1820}
1821
1823 size_t num_grids, void * grids_,
1824 void * buffer, int buffer_size, int * position, MPI_Comm comm) {
1825
1826 struct yac_couple_config_grid * grids = grids_;
1827
1828 YAC_ASSERT(
1829 num_grids <= INT_MAX,
1830 "ERROR(yac_couple_config_pack_grids): num_grids bigger than INT_MAX")
1831
1832 int num_grids_int = (int)num_grids;
1834 MPI_Pack(
1835 &num_grids_int, 1, MPI_INT, buffer, buffer_size, position, comm), comm);
1836
1837 for (size_t i = 0; i < num_grids; ++i)
1839 grids + i, buffer, buffer_size, position, comm);
1840}
1841
1843 struct yac_couple_config_field * field,
1844 void * buffer, int buffer_size, int * position, MPI_Comm comm) {
1845
1846 YAC_ASSERT(
1847 field->grid_idx <= INT_MAX,
1848 "ERROR(yac_couple_config_pack_field): grid_idx is too big")
1849 YAC_ASSERT(
1850 (field->collection_size < INT_MAX) ||
1851 (field->collection_size == SIZE_MAX),
1852 "ERROR(yac_couple_config_pack_field): invalid collection size")
1853
1854 int grid_idx = (int)field->grid_idx;
1855 double frac_mask_fallback_value = field->frac_mask_fallback_value;
1856 int collection_size =
1857 (field->collection_size == SIZE_MAX)?
1858 INT_MAX:(int)(field->collection_size);
1859
1861 MPI_Pack(
1862 &grid_idx, 1, MPI_INT, buffer, buffer_size, position, comm), comm);
1864 MPI_Pack(
1865 &frac_mask_fallback_value, 1, MPI_DOUBLE, buffer, buffer_size,
1866 position, comm), comm);
1868 MPI_Pack(
1869 &collection_size, 1, MPI_INT, buffer, buffer_size, position, comm), comm);
1870 yac_couple_config_pack_string(field->name, buffer,
1871 buffer_size, position, comm);
1873 buffer_size, position, comm);
1875 buffer_size, position, comm);
1876}
1877
1879 size_t num_fields, void * fields_, void * buffer, int buffer_size,
1880 int * position, MPI_Comm comm) {
1881
1882 struct yac_couple_config_field * fields = fields_;
1883
1884 YAC_ASSERT(
1885 num_fields <= INT_MAX,
1886 "ERROR(yac_couple_config_pack_fields): "
1887 "num_fields bigger than INT_MAX")
1888
1889 int num_fields_int = (int)num_fields;
1891 MPI_Pack(
1892 &num_fields_int, 1, MPI_INT,
1893 buffer, buffer_size, position, comm), comm);
1894
1895 for (size_t i = 0; i < num_fields; ++i)
1897 fields + i, buffer, buffer_size, position, comm);
1898}
1899
1901 struct yac_couple_config_component * component,
1902 void * buffer, int buffer_size, int * position, MPI_Comm comm) {
1903
1905 component->name, buffer, buffer_size, position, comm);
1906}
1907
1909 size_t num_components, void * components_,
1910 void * buffer, int buffer_size, int * position, MPI_Comm comm) {
1911
1912 struct yac_couple_config_component * components = components_;
1913
1914 YAC_ASSERT(
1915 num_components <= INT_MAX,
1916 "ERROR(yac_couple_config_pack_components): "
1917 "num_components bigger than INT_MAX")
1918
1919 int num_components_int = (int)num_components;
1921 MPI_Pack(
1922 &num_components_int, 1, MPI_INT, buffer, buffer_size, position, comm), comm);
1923
1924 for (size_t i = 0; i < num_components; ++i)
1926 components + i, buffer, buffer_size, position, comm);
1927}
1928
1930 struct yac_couple_config_field_couple * field_couple,
1931 void * buffer, int buffer_size, int * position, MPI_Comm comm) {
1932
1933 YAC_ASSERT(
1934 field_couple->source.component_idx <= INT_MAX,
1935 "ERROR(yac_couple_config_pack_field_couple): "
1936 "source.component_idx bigger than INT_MAX")
1937 YAC_ASSERT(
1938 field_couple->source.field_idx <= INT_MAX,
1939 "ERROR(yac_couple_config_pack_field_couple): "
1940 "source.field_idx bigger than INT_MAX")
1941 YAC_ASSERT(
1942 field_couple->target.component_idx <= INT_MAX,
1943 "ERROR(yac_couple_config_pack_field_couple): "
1944 "target.component_idx bigger than INT_MAX")
1945 YAC_ASSERT(
1946 field_couple->target.field_idx <= INT_MAX,
1947 "ERROR(yac_couple_config_pack_field_couple): "
1948 "target.field_idx bigger than INT_MAX")
1949 YAC_ASSERT(
1950 field_couple->mapping_on_source <= INT_MAX,
1951 "ERROR(yac_couple_config_pack_field_couple): "
1952 "mapping_on_source bigger than INT_MAX")
1953 YAC_ASSERT(
1954 field_couple->coupling_period_operation <= INT_MAX,
1955 "ERROR(yac_couple_config_pack_field_couple): "
1956 "coupling_period_operation bigger than INT_MAX")
1957 YAC_ASSERT(
1958 field_couple->enforce_write_weight_file <= INT_MAX,
1959 "ERROR(yac_couple_config_pack_field_couple): "
1960 "enforce_write_weight_file bigger than INT_MAX")
1961 YAC_ASSERT(
1962 field_couple->num_src_mask_names <= INT_MAX,
1963 "ERROR(yac_couple_config_pack_field_couple): "
1964 "num_src_mask_names bigger than INT_MAX")
1965
1966 int ints[10] = {
1967 field_couple->source.component_idx,
1968 field_couple->source.field_idx,
1969 field_couple->source.lag,
1970 field_couple->target.component_idx,
1971 field_couple->target.field_idx,
1972 field_couple->target.lag,
1973 field_couple->mapping_on_source,
1974 field_couple->coupling_period_operation,
1975 field_couple->enforce_write_weight_file,
1976 field_couple->num_src_mask_names};
1977
1979 MPI_Pack(ints, 10, MPI_INT, buffer, buffer_size, position, comm), comm);
1980
1982 field_couple->coupling_period, buffer, buffer_size, position, comm);
1984 field_couple->weight_file_name, buffer, buffer_size, position, comm);
1985
1986 double doubles[2] = {
1987 field_couple->scale_factor,
1988 field_couple->scale_summand};
1989
1991 MPI_Pack(doubles, 2, MPI_DOUBLE, buffer, buffer_size, position, comm),
1992 comm);
1993
1995 field_couple->interp_stack, buffer, buffer_size, position, comm);
1996
1997 if (field_couple->num_src_mask_names > 0)
1998 for (size_t i = 0; i < field_couple->num_src_mask_names; ++i)
2000 field_couple->src_mask_names[i], buffer, buffer_size, position, comm);
2001
2003 field_couple->tgt_mask_name, buffer, buffer_size, position, comm);
2004}
2005
2007 size_t num_field_couples, void * field_couples_,
2008 void * buffer, int buffer_size, int * position, MPI_Comm comm) {
2009
2010 struct yac_couple_config_field_couple * field_couples = field_couples_;
2011
2012 YAC_ASSERT(
2013 num_field_couples <= INT_MAX,
2014 "ERROR(yac_couple_config_pack_field_couples): "
2015 "num_field_couples bigger than INT_MAX")
2016
2017 int num_field_couples_int = (int)num_field_couples;
2019 MPI_Pack(
2020 &num_field_couples_int, 1, MPI_INT,
2021 buffer, buffer_size, position, comm), comm);
2022
2023 for (size_t i = 0; i < num_field_couples; ++i)
2025 field_couples + i, buffer, buffer_size, position, comm);
2026}
2027
2029 struct yac_couple_config_couple * couple,
2030 void * buffer, int buffer_size, int * position, MPI_Comm comm) {
2031
2032 YAC_ASSERT(
2033 (couple->component_indices[0] <= INT_MAX) &&
2034 (couple->component_indices[1] <= INT_MAX),
2035 "ERROR(yac_couple_config_pack_couple_basic): "
2036 "component_indices bigger than INT_MAX")
2037
2038 int component_indices[2] = {(int)(couple->component_indices[0]),
2039 (int)(couple->component_indices[1])};
2040
2042 MPI_Pack(
2043 component_indices, 2, MPI_INT, buffer, buffer_size, position, comm),
2044 comm);
2045}
2046
2048 size_t num_couples, void * couples_,
2049 void * buffer, int buffer_size, int * position, MPI_Comm comm) {
2050
2051 struct yac_couple_config_couple * couples = couples_;
2052
2053 YAC_ASSERT(
2054 num_couples <= INT_MAX,
2055 "ERROR(yac_couple_config_pack_couples_basic): "
2056 "num_couples bigger than INT_MAX")
2057
2058 int num_couples_int = (int)num_couples;
2060 MPI_Pack(
2061 &num_couples_int, 1, MPI_INT, buffer, buffer_size, position, comm), comm);
2062
2063 for (size_t i = 0; i < num_couples; ++i)
2065 couples + i, buffer, buffer_size, position, comm);
2066}
2067
2069 void * buffer, int buffer_size, int * position, MPI_Comm comm) {
2070
2071 int string_len;
2073 MPI_Unpack(
2074 buffer, buffer_size, position, &string_len, 1, MPI_INT, comm), comm);
2075
2076 if (string_len <= 0) return NULL;
2077
2078 char * string = xmalloc(((size_t)string_len + 1) * sizeof(*string));
2080 MPI_Unpack(
2081 buffer, buffer_size, position, string, string_len, MPI_CHAR, comm), comm);
2082 string[string_len] = '\0';
2083 return string;
2084}
2085
2087 void * buffer, int buffer_size, int * position,
2088 struct yac_couple_config_grid * grid, MPI_Comm comm) {
2089
2090 grid->name =
2091 yac_couple_config_unpack_string(buffer, buffer_size, position, comm);
2092 grid->metadata = NULL;
2093}
2094
2096 void * buffer, int buffer_size, int * position,
2097 size_t * num_grids, void * grids_, MPI_Comm comm) {
2098
2099 struct yac_couple_config_grid ** grids = grids_;
2100
2101 int num_grids_int;
2103 MPI_Unpack(
2104 buffer, buffer_size, position, &num_grids_int, 1, MPI_INT, comm), comm);
2105
2106 YAC_ASSERT(
2107 num_grids_int >= 0,
2108 "ERROR(yac_couple_config_unpack_grids): invalid number of grids")
2109
2110 *num_grids = (size_t)num_grids_int;
2111 *grids = xmalloc(*num_grids * sizeof(**grids));
2112
2113 for (size_t i = 0; i < *num_grids; ++i)
2115 buffer, buffer_size, position, *grids + i, comm);
2116}
2117
2119 void * buffer, int buffer_size, int * position,
2120 struct yac_couple_config_field * field, MPI_Comm comm) {
2121
2122 int grid_idx;
2123 double frac_mask_fallback_value;
2124 int collection_size;
2125
2127 MPI_Unpack(
2128 buffer, buffer_size, position, &grid_idx, 1, MPI_INT, comm), comm);
2130 MPI_Unpack(
2131 buffer, buffer_size, position, &frac_mask_fallback_value, 1,
2132 MPI_DOUBLE, comm), comm);
2134 MPI_Unpack(
2135 buffer, buffer_size, position, &collection_size, 1, MPI_INT, comm),
2136 comm);
2137
2138 YAC_ASSERT(
2139 grid_idx >= 0,
2140 "ERROR(yac_couple_config_unpack_field): invalid number of grid_idx")
2141 YAC_ASSERT(
2142 collection_size >= 0,
2143 "ERROR(yac_couple_config_unpack_field): invalid collection_size")
2144
2145 field->grid_idx = (size_t)grid_idx;
2146 field->frac_mask_fallback_value = frac_mask_fallback_value;
2147 field->collection_size =
2148 (collection_size == INT_MAX)?SIZE_MAX:((size_t)collection_size);
2149 field->name = yac_couple_config_unpack_string(buffer,
2150 buffer_size, position, comm);
2152 buffer_size, position, comm);
2154 buffer_size, position, comm);
2155}
2156
2158 void * buffer, int buffer_size, int * position,
2159 size_t * num_fields, void * fields_, MPI_Comm comm) {
2160
2161 struct yac_couple_config_field ** fields = fields_;
2162
2163 int num_fields_int;
2165 MPI_Unpack(
2166 buffer, buffer_size, position, &num_fields_int, 1,
2167 MPI_INT, comm), comm);
2168
2169 YAC_ASSERT(
2170 num_fields_int >= 0,
2171 "ERROR(yac_couple_config_unpack_fields): "
2172 "invalid number of fields")
2173
2174 *num_fields = (size_t)num_fields_int;
2175 *fields =
2176 xmalloc(*num_fields * sizeof(**fields));
2177
2178 for (size_t i = 0; i < *num_fields; ++i)
2180 buffer, buffer_size, position, *fields + i, comm);
2181}
2182
2184 void * buffer, int buffer_size, int * position,
2185 struct yac_couple_config_component * component, MPI_Comm comm) {
2186
2187 component->name =
2188 yac_couple_config_unpack_string(buffer, buffer_size, position, comm);
2189 component->metadata = NULL;
2190 component->num_fields = 0;
2191 component->fields = NULL;
2192}
2193
2195 void * buffer, int buffer_size, int * position,
2196 size_t * num_components, void * components_, MPI_Comm comm) {
2197
2198 struct yac_couple_config_component ** components = components_;
2199
2200 int num_components_int;
2202 MPI_Unpack(
2203 buffer, buffer_size, position, &num_components_int, 1, MPI_INT, comm),
2204 comm);
2205
2206 YAC_ASSERT(
2207 num_components_int >= 0,
2208 "ERROR(yac_couple_config_unpack_components): "
2209 "invalid number of components")
2210
2211 *num_components = (size_t)num_components_int;
2213
2214 for (size_t i = 0; i < *num_components; ++i)
2216 buffer, buffer_size, position, *components + i, comm);
2217}
2218
2220 void * buffer, int buffer_size, int * position,
2221 struct yac_couple_config_field_couple * field_couple,
2222 MPI_Comm comm) {
2223
2224 int ints[10];
2226 MPI_Unpack(
2227 buffer, buffer_size, position, ints, 10, MPI_INT, comm), comm);
2228
2229 YAC_ASSERT(
2230 ints[0] >= 0,
2231 "ERROR(yac_couple_config_unpack_field_couple): "
2232 "invalid source.component_idx")
2233 YAC_ASSERT(
2234 ints[1] >= 0,
2235 "ERROR(yac_couple_config_unpack_field_couple): "
2236 "invalid source.field_idx")
2237 YAC_ASSERT(
2238 ints[3] >= 0,
2239 "ERROR(yac_couple_config_unpack_field_couple): "
2240 "target.component_idx bigger than INT_MAX")
2241 YAC_ASSERT(
2242 ints[4] >= 0,
2243 "ERROR(yac_couple_config_unpack_field_couple): "
2244 "invalid target.field_idx")
2245 YAC_ASSERT(
2246 ints[6] >= 0,
2247 "ERROR(yac_couple_config_unpack_field_couple): "
2248 "invalid mapping_on_source")
2249 YAC_ASSERT(
2250 ints[7] >= 0,
2251 "ERROR(yac_couple_config_unpack_field_couple): "
2252 "invalid coupling_period_operation")
2253 YAC_ASSERT(
2254 ints[8] >= 0,
2255 "ERROR(yac_couple_config_unpack_field_couple): "
2256 "invalid enforce_write_weight_file")
2257 YAC_ASSERT(
2258 ints[9] >= 0,
2259 "ERROR(yac_couple_config_unpack_field_couple): "
2260 "invalid num_src_mask_names")
2261
2262 field_couple->source.component_idx = (size_t)(ints[0]);
2263 field_couple->source.field_idx = (size_t)(ints[1]);
2264 field_couple->source.lag = ints[2];
2265 field_couple->target.component_idx = (size_t)(ints[3]);
2266 field_couple->target.field_idx = (size_t)(ints[4]);
2267 field_couple->target.lag = ints[5];
2268 field_couple->mapping_on_source = ints[6];
2269 field_couple->coupling_period_operation =
2270 (enum yac_reduction_type)(ints[7]);
2271 field_couple->enforce_write_weight_file = ints[8];
2272 field_couple->num_src_mask_names = (size_t)(ints[9]);
2273
2274 field_couple->coupling_period =
2275 yac_couple_config_unpack_string(buffer, buffer_size, position, comm);
2276 field_couple->weight_file_name =
2277 yac_couple_config_unpack_string(buffer, buffer_size, position, comm);
2278
2279 double doubles[2];
2281 MPI_Unpack(
2282 buffer, buffer_size, position, doubles, 2, MPI_DOUBLE, comm), comm);
2283
2284 field_couple->scale_factor = doubles[0];
2285 field_couple->scale_summand = doubles[1];
2286
2287 field_couple->interp_stack =
2288 yac_interp_stack_config_unpack(buffer, buffer_size, position, comm);
2289
2290 if (field_couple->num_src_mask_names > 0) {
2291 field_couple->src_mask_names =
2292 xmalloc(
2293 field_couple->num_src_mask_names *
2294 sizeof(*(field_couple->src_mask_names)));
2295 for (size_t i = 0; i < field_couple->num_src_mask_names; ++i)
2296 field_couple->src_mask_names[i] =
2298 buffer, buffer_size, position, comm);
2299 } else {
2300 field_couple->src_mask_names = NULL;
2301 }
2302
2303 field_couple->tgt_mask_name =
2305 buffer, buffer_size, position, comm);
2306}
2307
2309 void * buffer, int buffer_size, int * position,
2310 size_t * num_field_couples, void * field_couples_, MPI_Comm comm) {
2311
2312 struct yac_couple_config_field_couple ** field_couples = field_couples_;
2313
2314 int num_field_couples_int;
2316 MPI_Unpack(
2317 buffer, buffer_size, position, &num_field_couples_int, 1,
2318 MPI_INT, comm), comm);
2319
2320 YAC_ASSERT(
2321 num_field_couples_int >= 0,
2322 "ERROR(yac_couple_config_unpack_field_couples): "
2323 "invalid number of field_couples")
2324
2325 *num_field_couples = (size_t)num_field_couples_int;
2326 *field_couples =
2327 xmalloc(*num_field_couples * sizeof(**field_couples));
2328
2329 for (size_t i = 0; i < *num_field_couples; ++i)
2331 buffer, buffer_size, position, *field_couples + i, comm);
2332}
2333
2335 void * buffer, int buffer_size, int * position,
2336 struct yac_couple_config_couple * couple, MPI_Comm comm) {
2337
2338 int component_indices[2];
2340 MPI_Unpack(
2341 buffer, buffer_size, position, component_indices, 2, MPI_INT, comm),
2342 comm);
2343
2344 YAC_ASSERT(
2345 (component_indices[0] >= 0) && (component_indices[1] >= 0),
2346 "ERROR(yac_couple_config_unpack_couple_basic): invalid component indices")
2347
2348 couple->component_indices[0] = (size_t)(component_indices[0]);
2349 couple->component_indices[1] = (size_t)(component_indices[1]);
2350 couple->num_field_couples = 0;
2351 couple->field_couples = NULL;
2352}
2353
2355 void * buffer, int buffer_size, int * position,
2356 size_t * num_couples, void * couples_, MPI_Comm comm) {
2357
2358 struct yac_couple_config_couple ** couples = couples_;
2359
2360 int num_couples_int;
2362 MPI_Unpack(
2363 buffer, buffer_size, position, &num_couples_int, 1, MPI_INT, comm),
2364 comm);
2365
2366 YAC_ASSERT(
2367 num_couples_int >= 0,
2368 "ERROR(yac_couple_config_unpack_couples_basic): "
2369 "invalid number of couples")
2370
2371 *num_couples = (size_t)num_couples_int;
2372 *couples = xmalloc(*num_couples * sizeof(**couples));
2373
2374 for (size_t i = 0; i < *num_couples; ++i)
2376 buffer, buffer_size, position, *couples + i, comm);
2377}
2378
2380 struct yac_couple_config * couple_config,
2381 char const * src_comp_name, char const * src_grid_name, char const * src_field_name,
2382 char const * tgt_comp_name, char const * tgt_grid_name, char const * tgt_field_name,
2383 char const * coupling_period, int time_reduction,
2384 struct yac_interp_stack_config * interp_stack,
2385 int src_lag, int tgt_lag,
2386 const char* weight_file_name, int mapping_on_source,
2387 double scale_factor, double scale_summand,
2388 size_t num_src_mask_names, char const * const * src_mask_names,
2389 char const * tgt_mask_name) {
2390
2391 YAC_ASSERT(src_comp_name && src_comp_name[0] != '\0',
2392 "ERROR(yac_couple_config_def_couple): invalid parameter: src_comp_name");
2393 YAC_ASSERT(src_grid_name && src_grid_name[0] != '\0',
2394 "ERROR(yac_couple_config_def_couple): invalid parameter: src_grid_name");
2395 YAC_ASSERT(src_field_name && src_field_name[0] != '\0',
2396 "ERROR(yac_couple_config_def_couple): invalid parameter: src_field_name");
2397 YAC_ASSERT(tgt_comp_name && tgt_comp_name[0] != '\0',
2398 "ERROR(yac_couple_config_def_couple): invalid parameter: tgt_comp_name");
2399 YAC_ASSERT(tgt_grid_name && tgt_grid_name[0] != '\0',
2400 "ERROR(yac_couple_config_def_couple): invalid parameter: tgt_grid_name");
2401 YAC_ASSERT(tgt_field_name && tgt_field_name[0] != '\0',
2402 "ERROR(yac_couple_config_def_couple): invalid parameter: tgt_field_name");
2403 YAC_ASSERT(coupling_period && coupling_period[0] != '\0',
2404 "ERROR(yac_couple_config_def_couple): invalid parameter: coupling_period");
2405 YAC_ASSERT(
2406 (time_reduction == TIME_NONE) ||
2407 (time_reduction == TIME_ACCUMULATE) ||
2408 (time_reduction == TIME_AVERAGE) ||
2409 (time_reduction == TIME_MINIMUM) ||
2410 (time_reduction == TIME_MAXIMUM),
2411 "ERROR(yac_couple_config_def_couple): invalid parameter: time_reduction");
2412 YAC_ASSERT_F(isnormal(scale_factor),
2413 "ERROR(yac_couple_config_def_couple): \"%lf\" is not a valid scale factor",
2414 scale_factor);
2415 YAC_ASSERT_F(isnormal(scale_summand) || (scale_summand == 0.0),
2416 "ERROR(yac_couple_config_def_couple): \"%lf\" is not a valid scale summand",
2417 scale_summand);
2418
2419 // get component indices
2420 size_t src_comp_idx =
2421 yac_couple_config_add_component_(couple_config, src_comp_name);
2422 size_t tgt_comp_idx =
2423 yac_couple_config_add_component_(couple_config, tgt_comp_name);
2424 size_t src_grid_idx =
2425 yac_couple_config_add_grid_(couple_config, src_grid_name);
2426 size_t tgt_grid_idx =
2427 yac_couple_config_add_grid_(couple_config, tgt_grid_name);
2428
2429 // check if couple exists
2430 size_t component_indices[2];
2431 if(src_comp_idx < tgt_comp_idx){
2432 component_indices[0] = src_comp_idx;
2433 component_indices[1] = tgt_comp_idx;
2434 }else{
2435 component_indices[0] = tgt_comp_idx;
2436 component_indices[1] = src_comp_idx;
2437 }
2438 struct yac_couple_config_couple * couple = NULL;
2439 for(size_t i = 0; (i < couple_config->num_couples) && !couple; ++i)
2440 if(couple_config->couples[i].component_indices[0] == component_indices[0] &&
2441 couple_config->couples[i].component_indices[1] == component_indices[1])
2442 couple = &couple_config->couples[i];
2443
2444 // create if couple does not exists
2445 if(!couple){
2446 couple_config->couples =
2447 xrealloc(
2448 couple_config->couples,
2449 (couple_config->num_couples + 1) * sizeof(*couple_config->couples));
2450 couple = &couple_config->couples[couple_config->num_couples];
2451 couple_config->num_couples++;
2452 couple->component_indices[0] = component_indices[0];
2453 couple->component_indices[1] = component_indices[1];
2454 couple->num_field_couples = 0;
2455 couple->field_couples = NULL;
2456 }
2457
2458 // get field indices
2459 size_t src_field_idx =
2461 couple_config, src_comp_idx, src_grid_idx,
2462 src_field_name, NULL, SIZE_MAX);
2463 size_t tgt_field_idx =
2465 couple_config, tgt_comp_idx, tgt_grid_idx,
2466 tgt_field_name, NULL, SIZE_MAX);
2467
2468 // check if field_couple exists
2469 for(size_t i = 0; i < couple->num_field_couples; ++i)
2471 !(couple->field_couples[i].source.component_idx == src_comp_idx &&
2472 couple->field_couples[i].target.component_idx == tgt_comp_idx &&
2473 couple->field_couples[i].source.field_idx == src_field_idx &&
2474 couple->field_couples[i].target.field_idx == tgt_field_idx),
2475 "ERROR(yac_couple_config_def_couple): Coupling is already defined \n"
2476 "source (comp_name: \"%s\" grid_name: \"%s\" field_name: \"%s\")\n"
2477 "target (comp_name: \"%s\" grid_name: \"%s\" field_name: \"%s\")\n",
2478 src_comp_name, src_grid_name, src_field_name,
2479 tgt_comp_name, tgt_grid_name, tgt_field_name);
2480
2481 couple->field_couples = xrealloc(couple->field_couples,
2482 (couple->num_field_couples + 1) * sizeof(*couple->field_couples));
2483 struct yac_couple_config_field_couple * field_couple =
2484 couple->field_couples + couple->num_field_couples;
2485 couple->num_field_couples++;
2486 field_couple->source.component_idx = src_comp_idx;
2487 field_couple->source.field_idx = src_field_idx;
2488 field_couple->source.lag = src_lag;
2489 field_couple->target.component_idx = tgt_comp_idx;
2490 field_couple->target.field_idx = tgt_field_idx;
2491 field_couple->target.lag = tgt_lag;
2492 field_couple->coupling_period = strdup(coupling_period);
2493 field_couple->coupling_period_operation =
2494 (enum yac_reduction_type)time_reduction;
2496 field_couple->mapping_on_source = mapping_on_source;
2497 field_couple->weight_file_name =
2499 field_couple->scale_factor = scale_factor;
2500 field_couple->scale_summand = scale_summand;
2501 field_couple->enforce_write_weight_file = weight_file_name != NULL;
2502 field_couple->num_src_mask_names = num_src_mask_names;
2503 if (num_src_mask_names > 0) {
2504 field_couple->src_mask_names =
2506 for (size_t i = 0; i < num_src_mask_names; ++i)
2507 field_couple->src_mask_names[i] = strdup(src_mask_names[i]);
2508 } else {
2509 field_couple->src_mask_names = NULL;
2510 }
2511 field_couple->tgt_mask_name =
2512 (tgt_mask_name != NULL)?strdup(tgt_mask_name):NULL;
2513}
2514
2516 struct yac_couple_config * couple_config, MPI_Comm comm) {
2517
2520 "start_datetime", (char**)&(couple_config->start_datetime), comm);
2522 "end_datetime", (char**)&(couple_config->end_datetime), comm);
2523}
2524
2526 struct yac_couple_config * couple_config, MPI_Comm comm) {
2527
2529 "config_filename", (char**)&(couple_config->config_filename), comm);
2530 couple_config_sync_config_filetype(&(couple_config->config_filetype), comm);
2531}
2532
2534 size_t(*get_pack_size)(size_t, void*, MPI_Comm);
2535 void(*pack)(size_t, void*, void *, int, int *, MPI_Comm);
2536 void(*unpack)(void *, int, int *, size_t *, void *, MPI_Comm);
2537 int(*compare)(void const*, void const*);
2538 void(*merge)(void*, void *, MPI_Comm);
2539 void(*free_data)(void*);
2575
2576static void dist_merge(size_t* len, void** arr, size_t element_size, MPI_Comm comm,
2577 struct dist_merge_vtable* vtable, size_t** idx_old_to_new){
2578 int rank;
2579 MPI_Comm_rank(comm, &rank);
2580
2581 // initialise
2582 unsigned char * input = *arr;
2583 size_t input_len = *len;
2584 *idx_old_to_new = xmalloc(input_len * sizeof(size_t));
2585 size_t* idx = xmalloc(input_len * sizeof(size_t));
2586 for(size_t i = 0;i<input_len;++i) idx[i] = i;
2587 unsigned char * arr_new = NULL;
2588 size_t len_new = 0;
2589
2590 // sort input data
2591 yac_qsort_index(input, input_len, element_size, vtable->compare, idx);
2592
2593 void * buffer = NULL;
2594 while(1){
2595
2596 // determine pack size of local data
2597 size_t pack_size =
2598 (input_len > 0)?vtable->get_pack_size(input_len, input, comm):0;
2599 YAC_ASSERT(
2600 pack_size <= LONG_MAX, "ERROR(dist_merge): packing size too big");
2601
2602 // determine rank with most amount of data
2603 struct {
2604 long pack_size;
2605 int rank;
2606 } data_pair = {.pack_size = (long)pack_size, .rank = rank};
2608 MPI_Allreduce(
2609 MPI_IN_PLACE, &data_pair, 1, MPI_LONG_INT, MPI_MAXLOC, comm), comm);
2610
2611 // if there is no more data to exchange
2612 if (data_pair.pack_size == 0) break;
2613
2614 // pack data into buffer buffer
2615 pack_size = (size_t)data_pair.pack_size;
2616 if (!buffer) buffer = xmalloc(pack_size);
2617 int position = 0;
2618 if(data_pair.rank == rank)
2619 vtable->pack(input_len, input, buffer, pack_size, &position, comm);
2620
2621 // broadcast and unpack data
2623 MPI_Bcast(buffer, pack_size, MPI_PACKED, data_pair.rank, comm), comm);
2624 void * recved = NULL;
2625 size_t num_recved;
2626 position = 0;
2627 vtable->unpack(buffer, pack_size, &position, &num_recved, &recved, comm);
2628
2629 // add received elements to list
2630 arr_new = xrealloc(arr_new, (len_new + num_recved)*element_size);
2631 memcpy(
2632 arr_new + len_new*element_size, recved, num_recved*element_size);
2633 free(recved);
2634
2635 // for all received elements
2636 size_t input_idx, input_len_new, i = len_new;
2637 len_new += num_recved;
2638 for(input_idx = 0, input_len_new = 0; i < len_new; ++i) {
2639 void* recved_element = arr_new + i*element_size;
2640
2641 // search for matching element in input list
2642 int cmp = 0;
2643 void * input_element = NULL;
2644 while ((input_idx < input_len) &&
2645 (((cmp =
2646 vtable->compare(
2647 ((input_element = input + input_idx * element_size)),
2648 recved_element))) < 0)) {
2649
2650 if (input_idx != input_len_new) {
2651 memcpy(
2652 input + input_len_new * element_size,
2653 input_element, element_size);
2654 idx[input_len_new] = idx[input_idx];
2655 }
2656 ++input_len_new;
2657 ++input_idx;
2658 }
2659
2660 // if the end of the local input list was reached
2661 if (input_idx == input_len) break;
2662
2663 // if a matching element was found in the input list
2664 if (!cmp) {
2665
2666 // merge input list element with received element
2667 if (vtable->merge)
2668 vtable->merge(recved_element, input_element, comm);
2669 // remove element
2670 vtable->free_data(input_element);
2671 // update index
2672 (*idx_old_to_new)[idx[input_idx]] = i;
2673 // upate input idx
2674 input_idx++;
2675
2676 // if no matching element was found in the input list
2677 } else if (vtable->merge) {
2678 // since the merge operation can potentially be collective, we have
2679 // to call it, even if no matching element was found
2680 vtable->merge(recved_element, NULL, comm);
2681 }
2682 }
2683 // process remaining received elements
2684 if (vtable->merge)
2685 for(; i < len_new; ++i)
2686 vtable->merge(arr_new + i*element_size, NULL, comm);
2687 // pack remaining elements in input list
2688 for (; input_idx < input_len; ++input_idx, ++input_len_new) {
2689 if (input_idx != input_len_new) {
2690 memcpy(
2691 input + input_len_new * element_size,
2692 input + input_idx * element_size, element_size);
2693 idx[input_len_new] = idx[input_idx];
2694 }
2695 }
2696 // update length of input list
2697 input_len = input_len_new;
2698 }
2699 free(buffer);
2700 free(input);
2701 free(idx);
2702 *arr = arr_new;
2703 *len = len_new;
2704}
2705
2706static void merge_grids(
2707 struct yac_couple_config * couple_config, MPI_Comm comm) {
2708
2709 size_t* old_to_new_idx;
2710 void * p_grids = couple_config->grids;
2711 dist_merge(
2712 &couple_config->num_grids, &p_grids,
2713 sizeof(couple_config->grids[0]),
2714 comm, &dist_merge_vtable_grid, &old_to_new_idx);
2715 couple_config->grids = p_grids;
2716
2717 // set new grid_idx in fields
2718 for(size_t comp_idx = 0; comp_idx < couple_config->num_components;
2719 ++comp_idx) {
2720 struct yac_couple_config_component * component =
2721 couple_config->components + comp_idx;
2722 for(size_t field_idx = 0; field_idx < component->num_fields; ++field_idx)
2723 component->fields[field_idx].grid_idx =
2724 old_to_new_idx[component->fields[field_idx].grid_idx];
2725 }
2726
2727 free(old_to_new_idx);
2728}
2729
2730static void merge_fields(
2731 struct yac_couple_config * couple_config, size_t comp_idx, MPI_Comm comm) {
2732
2733 struct yac_couple_config_component * component =
2734 couple_config->components + comp_idx;
2735
2736 size_t* old_to_new_idx;
2737 void * p_fields = component->fields;
2738 dist_merge(
2739 &component->num_fields, &p_fields,
2740 sizeof(component->fields[0]),
2741 comm, &dist_merge_vtable_field, &old_to_new_idx);
2742 component->fields = p_fields;
2743
2744 // set new field_idx in all field_couples
2745 for(size_t couple_idx = 0; couple_idx < couple_config->num_couples;
2746 ++couple_idx) {
2747 struct yac_couple_config_couple * couple =
2748 couple_config->couples + couple_idx;
2749 for(size_t field_couple_idx = 0;
2750 field_couple_idx < couple->num_field_couples; ++field_couple_idx){
2751 struct yac_couple_config_field_couple * field_couple =
2752 couple->field_couples + field_couple_idx;
2753 if(field_couple->source.component_idx == comp_idx)
2754 field_couple->source.field_idx =
2755 old_to_new_idx[field_couple->source.field_idx];
2756 if(field_couple->target.component_idx == comp_idx)
2757 field_couple->target.field_idx =
2758 old_to_new_idx[field_couple->target.field_idx];
2759 }
2760 }
2761
2762 free(old_to_new_idx);
2763}
2764
2766 struct yac_couple_config * couple_config, MPI_Comm comm) {
2767
2768 // distribute and merge basic component information while keeping the
2769 // individual fields
2770 size_t* old_to_new_idx;
2771 void * p_components = couple_config->components;
2772 dist_merge(
2773 &couple_config->num_components, &p_components,
2774 sizeof(couple_config->components[0]),
2775 comm, &dist_merge_vtable_component, &old_to_new_idx);
2776 couple_config->components = p_components;
2777
2778 // set new component_idx in couples
2779 for(size_t couple_idx = 0; couple_idx < couple_config->num_couples;
2780 ++couple_idx) {
2781 struct yac_couple_config_couple * couple =
2782 couple_config->couples + couple_idx;
2783 couple->component_indices[0] = old_to_new_idx[couple->component_indices[0]];
2784 couple->component_indices[1] = old_to_new_idx[couple->component_indices[1]];
2785 for(size_t field_couple_idx = 0;
2786 field_couple_idx < couple->num_field_couples; ++field_couple_idx) {
2787 couple->field_couples[field_couple_idx].source.component_idx =
2788 old_to_new_idx[
2789 couple->field_couples[field_couple_idx].source.component_idx];
2790 couple->field_couples[field_couple_idx].target.component_idx =
2791 old_to_new_idx[
2792 couple->field_couples[field_couple_idx].target.component_idx];
2793 }
2794 }
2795 free(old_to_new_idx);
2796
2797 // merge the fields of each component
2798 for (size_t comp_idx = 0; comp_idx < couple_config->num_components; ++comp_idx)
2799 merge_fields(couple_config, comp_idx, comm);
2800}
2801
2803 size_t * num_field_couples,
2804 struct yac_couple_config_field_couple ** field_couples, MPI_Comm comm) {
2805
2806 size_t* old_to_new_idx;
2807 void * p_field_couples = *field_couples;
2808 dist_merge(
2809 num_field_couples, &p_field_couples, sizeof(**field_couples),
2810 comm, &dist_merge_vtable_field_couple, &old_to_new_idx);
2811 free(old_to_new_idx);
2812 *field_couples = p_field_couples;
2813}
2814
2815static void merge_couples(
2816 struct yac_couple_config * couple_config, MPI_Comm comm) {
2817
2818 // distribute and merge basic couple information while keeping the
2819 // individual field couples
2820 size_t* old_to_new_idx;
2821 void * p_couples = couple_config->couples;
2822 dist_merge(
2823 &couple_config->num_couples, &p_couples,
2824 sizeof(couple_config->couples[0]),
2825 comm, &dist_merge_vtable_couple, &old_to_new_idx);
2826 couple_config->couples = p_couples;
2827 free(old_to_new_idx);
2828}
2829
2831 struct yac_couple_config * couple_config, MPI_Comm comm){
2832
2833 // sync time stuff
2834 couple_config_sync_time(couple_config, comm);
2835
2836 // sync output configuration
2837 couple_config_sync_config_file(couple_config, comm);
2838
2839 merge_grids(couple_config, comm);
2840 merge_components(couple_config, comm);
2841 merge_couples(couple_config, comm);
2842
2843 if (couple_config->config_filetype != YAC_CONFIG_FILETYPE_INVALID) {
2844
2845 int rank;
2846 MPI_Comm_rank(comm, &rank);
2847
2848 // rank 0 writes the coupling configuration to file
2849 if (rank == 0) {
2850
2851 FILE * config_file = fopen(couple_config->config_filename, "w");
2852
2854 config_file != NULL,
2855 "ERROR(yac_couple_config_sync): "
2856 "failed to create coupling configuration file \"%s\"",
2857 couple_config->config_filename);
2858
2860 (couple_config->config_filetype == YAC_CONFIG_FILETYPE_YAML) ||
2861 (couple_config->config_filetype == YAC_CONFIG_FILETYPE_JSON),
2862 "ERROR(yac_couple_config_sync): "
2863 "invalid coupling configuration filetype (type = %d)",
2864 couple_config->config_filetype);
2865
2866 int emit_flags;
2867 switch(couple_config->config_filetype) {
2868 default:
2870 emit_flags = YAC_YAML_EMITTER_DEFAULT;
2871 break;
2873 emit_flags = YAC_YAML_EMITTER_JSON;
2874 break;
2875 }
2876
2877 // char * str_couple_config =
2878 // yac_yaml_emit_coupling(couple_config, emit_flags);
2879 char * str_couple_config = NULL;
2880 str_couple_config =
2881 yac_yaml_emit_coupling(couple_config, emit_flags);
2882
2883 fputs(str_couple_config, config_file);
2884 free(str_couple_config);
2885 fclose(config_file);
2886 }
2887 yac_mpi_call(MPI_Barrier(comm), comm);
2888 }
2889}
2890
2892 struct yac_couple_config * couple_config, char const * filename,
2893 enum yac_config_filetype filetype) {
2894
2895 YAC_ASSERT(
2896 filename != NULL,
2897 "ERROR(yac_couple_config_enable_output): filename is NULL");
2898
2900 (filetype == YAC_CONFIG_FILETYPE_YAML) ||
2901 (filetype == YAC_CONFIG_FILETYPE_JSON) ||
2902 (filetype == YAC_CONFIG_FILETYPE_INVALID),
2903 "ERROR(yac_couple_config_enable_output): "
2904 "unsupported output configuration filetype (type = %d)",
2905 (int)filetype);
2906
2908 (filetype == YAC_CONFIG_FILETYPE_YAML) ||
2909 (filetype == YAC_CONFIG_FILETYPE_JSON),
2910 "ERROR(yac_couple_config_enable_output): "
2911 "invalid output configuration filetype (type = %d)",
2912 (int)filetype);
2913
2914 free((void*)(couple_config->config_filename));
2915
2916 couple_config->config_filename = strdup(filename);
2917 couple_config->config_filetype = filetype;
2918}
2919
2921 struct yac_couple_config * couple_config) {
2922
2923 free((void*)(couple_config->config_filename));
2924 couple_config->config_filename = NULL;
2926}
int YAC_YAML_EMITTER_JSON
emit to JSON format
Definition config_yaml.c:62
int YAC_YAML_EMITTER_DEFAULT
emit to YAML format
Definition config_yaml.c:61
char * yac_yaml_emit_coupling(struct yac_couple_config *couple_config, int emit_flags)
#define UNUSED(x)
Definition core.h:71
void yac_couple_config_grid_set_metadata(struct yac_couple_config *couple_config, char const *grid_name, const char *metadata)
static size_t yac_couple_config_get_field_pack_size(struct yac_couple_config_field *field, MPI_Comm comm)
static void check_field_couple_idx(struct yac_couple_config *couple_config, size_t couple_idx, size_t field_couple_idx, char const *routine_name, int line)
static void yac_couple_config_unpack_field_couples(void *buffer, int buffer_size, int *position, size_t *num_field_couples, void *field_couples_, MPI_Comm comm)
static void yac_couple_config_grid_merge(void *a_, void *b_, MPI_Comm comm)
int yac_couple_config_field_is_valid(struct yac_couple_config *couple_config, size_t component_idx, size_t field_idx)
static void dist_merge(size_t *len, void **arr, size_t element_size, MPI_Comm comm, struct dist_merge_vtable *vtable, size_t **idx_old_to_new)
static void couple_config_sync_time(struct yac_couple_config *couple_config, MPI_Comm comm)
static void yac_couple_config_field_couple_free(void *field_couple_)
static int yac_couple_config_grid_compare(void const *a, void const *b)
static size_t yac_couple_config_get_grid_pack_size(void *grid_, MPI_Comm comm)
static void yac_couple_config_component_merge(void *a_, void *b_, MPI_Comm comm)
void yac_couple_config_field_enable_frac_mask(struct yac_couple_config *couple_config, char const *comp_name, char const *grid_name, char const *field_name, double frac_mask_fallback_value)
static void yac_couple_config_unpack_component(void *buffer, int buffer_size, int *position, struct yac_couple_config_component *component, MPI_Comm comm)
static void couple_config_sync_string(char const *string_name, char **string, MPI_Comm comm)
static void merge_fields(struct yac_couple_config *couple_config, size_t comp_idx, MPI_Comm comm)
static char const * string_dup(char const *string)
static void yac_couple_config_couple_free(void *couple_)
static void yac_couple_config_pack_components(size_t num_components, void *components_, void *buffer, int buffer_size, int *position, MPI_Comm comm)
size_t yac_couple_config_get_num_fields(struct yac_couple_config *couple_config, size_t component_idx)
void yac_couple_config_get_field_grid_names(struct yac_couple_config *couple_config, size_t couple_idx, size_t field_couple_idx, char const **src_grid_name, char const **tgt_grid_name)
int yac_couple_config_contains_grid_name(struct yac_couple_config *couple_config, char const *grid_name)
static void yac_couple_config_pack_field_couples(size_t num_field_couples, void *field_couples_, void *buffer, int buffer_size, int *position, MPI_Comm comm)
void yac_couple_config_sync(struct yac_couple_config *couple_config, MPI_Comm comm)
void yac_couple_config_set_datetime(struct yac_couple_config *couple_config, char const *start, char const *end)
static void yac_couple_config_unpack_couple_basic(void *buffer, int buffer_size, int *position, struct yac_couple_config_couple *couple, MPI_Comm comm)
void yac_couple_config_component_add_field(struct yac_couple_config *couple_config, const char *component_name, const char *grid_name, const char *name, char const *timestep, size_t collection_size)
char const * yac_couple_config_get_field_timestep(struct yac_couple_config *couple_config, char const *component_name, char const *grid_name, char const *field_name)
static void yac_couple_config_pack_couples_basic(size_t num_couples, void *couples_, void *buffer, int buffer_size, int *position, MPI_Comm comm)
void yac_couple_config_disable_output(struct yac_couple_config *couple_config)
int yac_couple_config_mapping_on_source(struct yac_couple_config *couple_config, size_t couple_idx, size_t field_couple_idx)
int yac_couple_config_get_field_role(struct yac_couple_config *couple_config, char const *component_name, char const *grid_name, char const *field_name)
int yac_couple_config_component_name_is_valid(struct yac_couple_config *couple_config, char const *component_name)
static void yac_couple_config_pack_couple_basic(struct yac_couple_config_couple *couple, void *buffer, int buffer_size, int *position, MPI_Comm comm)
static size_t yac_couple_config_get_field_couples_pack_size(size_t num_field_couples, void *field_couples_, MPI_Comm comm)
static char * yac_couple_config_unpack_string(void *buffer, int buffer_size, int *position, MPI_Comm comm)
char const * yac_couple_config_get_coupling_period(struct yac_couple_config *couple_config, size_t couple_idx, size_t field_couple_idx)
struct dist_merge_vtable dist_merge_vtable_component
void yac_couple_config_field_set_metadata(struct yac_couple_config *couple_config, const char *comp_name, const char *grid_name, const char *field_name, const char *metadata)
struct yac_interp_stack_config * yac_couple_config_get_interp_stack(struct yac_couple_config *couple_config, size_t couple_idx, size_t field_couple_idx)
static void yac_couple_config_unpack_grid(void *buffer, int buffer_size, int *position, struct yac_couple_config_grid *grid, MPI_Comm comm)
static size_t yac_couple_config_get_field_couple_pack_size(struct yac_couple_config_field_couple *field_couple, MPI_Comm comm)
char const * yac_couple_config_get_field_name(struct yac_couple_config *couple_config, size_t component_idx, size_t field_idx)
size_t yac_couple_config_get_field_idx(struct yac_couple_config *couple_config, size_t component_idx, size_t grid_idx, char const *field_name)
static size_t yac_couple_config_get_components_pack_size(size_t num_components, void *components_, MPI_Comm comm)
static void couple_config_sync_calendar(MPI_Comm comm)
static size_t yac_couple_config_add_component_(struct yac_couple_config *couple_config, char const *name)
void yac_couple_config_component_set_metadata(struct yac_couple_config *couple_config, char const *comp_name, const char *metadata)
static size_t yac_couple_config_add_grid_(struct yac_couple_config *couple_config, char const *name)
struct dist_merge_vtable dist_merge_vtable_field
static void merge_components(struct yac_couple_config *couple_config, MPI_Comm comm)
static void check_couple_idx(struct yac_couple_config *couple_config, size_t couple_idx, char const *routine_name, int line)
size_t yac_couple_config_get_num_couple_fields(struct yac_couple_config *couple_config, size_t couple_idx)
void yac_couple_config_get_src_mask_names(struct yac_couple_config *couple_config, size_t couple_idx, size_t field_couple_idx, char const *const **mask_names, size_t *num_mask_names)
static int yac_couple_config_component_compare(void const *a_, void const *b_)
static void yac_couple_config_unpack_field(void *buffer, int buffer_size, int *position, struct yac_couple_config_field *field, MPI_Comm comm)
static void yac_couple_config_unpack_couples_basic(void *buffer, int buffer_size, int *position, size_t *num_couples, void *couples_, MPI_Comm comm)
void yac_couple_config_enable_output(struct yac_couple_config *couple_config, char const *filename, enum yac_config_filetype filetype)
size_t yac_couple_config_get_num_couples(struct yac_couple_config *couple_config)
static size_t yac_couple_config_get_couple_pack_size_basic(struct yac_couple_config_couple *couple, MPI_Comm comm)
static int yac_couple_config_field_compare(void const *a_, void const *b_)
static void yac_couple_config_couple_merge(void *a_, void *b_, MPI_Comm comm)
int yac_couple_config_get_target_lag(struct yac_couple_config *couple_config, size_t couple_idx, size_t field_couple_idx)
static void yac_couple_config_pack_field_couple(struct yac_couple_config_field_couple *field_couple, void *buffer, int buffer_size, int *position, MPI_Comm comm)
double yac_couple_config_get_frac_mask_fallback_value(struct yac_couple_config *couple_config, char const *component_name, char const *grid_name, char const *field_name)
static size_t yac_couple_config_get_fields_pack_size(size_t num_fields, void *fields_, MPI_Comm comm)
int yac_couple_config_get_source_lag(struct yac_couple_config *couple_config, size_t couple_idx, size_t field_couple_idx)
static int yac_couple_config_couple_compare_basic(void const *a_, void const *b_)
static void couple_config_sync_config_filetype(enum yac_config_filetype *filetype, MPI_Comm comm)
const char * yac_couple_config_field_get_metadata(struct yac_couple_config *couple_config, const char *comp_name, const char *grid_name, const char *field_name)
static void merge_grids(struct yac_couple_config *couple_config, MPI_Comm comm)
void yac_couple_config_get_couple_component_names(struct yac_couple_config *couple_config, size_t couple_idx, char const *couple_component_names[2])
const char * yac_couple_config_grid_get_metadata(struct yac_couple_config *couple_config, const char *grid_name)
static void merge_field_couples(size_t *num_field_couples, struct yac_couple_config_field_couple **field_couples, MPI_Comm comm)
static void yac_couple_config_unpack_grids(void *buffer, int buffer_size, int *position, size_t *num_grids, void *grids_, MPI_Comm comm)
const char * yac_couple_config_component_get_metadata(struct yac_couple_config *couple_config, const char *comp_name)
double yac_couple_config_get_scale_factor(struct yac_couple_config *couple_config, size_t couple_idx, size_t field_couple_idx)
char const * yac_couple_config_get_end_datetime(struct yac_couple_config *couple_config)
static void yac_couple_config_pack_fields(size_t num_fields, void *fields_, void *buffer, int buffer_size, int *position, MPI_Comm comm)
size_t yac_couple_config_get_grid_idx(struct yac_couple_config *couple_config, char const *grid_name)
size_t yac_couple_config_get_component_idx(struct yac_couple_config *couple_config, char const *component_name)
char const * yac_couple_config_get_tgt_mask_name(struct yac_couple_config *couple_config, size_t couple_idx, size_t field_couple_idx)
void yac_couple_config_add_grid(struct yac_couple_config *couple_config, char const *name)
char const * yac_couple_config_get_grid_name(struct yac_couple_config *couple_config, size_t grid_idx)
double yac_couple_config_get_scale_summand(struct yac_couple_config *couple_config, size_t couple_idx, size_t field_couple_idx)
static void yac_couple_config_unpack_fields(void *buffer, int buffer_size, int *position, size_t *num_fields, void *fields_, MPI_Comm comm)
enum yac_reduction_type yac_couple_config_get_coupling_period_operation(struct yac_couple_config *couple_config, size_t couple_idx, size_t field_couple_idx)
static void yac_couple_config_pack_grids(size_t num_grids, void *grids_, void *buffer, int buffer_size, int *position, MPI_Comm comm)
static size_t yac_couple_config_get_string_pack_size(char const *string, MPI_Comm comm)
int yac_couple_config_enforce_write_weight_file(struct yac_couple_config *couple_config, size_t couple_idx, size_t field_couple_idx)
static void yac_couple_config_grid_free(void *grid_)
static size_t yac_couple_config_get_couples_pack_size_basic(size_t num_couples, void *couples_, MPI_Comm comm)
static void yac_couple_config_unpack_field_couple(void *buffer, int buffer_size, int *position, struct yac_couple_config_field_couple *field_couple, MPI_Comm comm)
void yac_couple_config_def_couple(struct yac_couple_config *couple_config, char const *src_comp_name, char const *src_grid_name, char const *src_field_name, char const *tgt_comp_name, char const *tgt_grid_name, char const *tgt_field_name, char const *coupling_period, int time_reduction, struct yac_interp_stack_config *interp_stack, int src_lag, int tgt_lag, const char *weight_file_name, int mapping_on_source, double scale_factor, double scale_summand, size_t num_src_mask_names, char const *const *src_mask_names, char const *tgt_mask_name)
static size_t yac_couple_config_get_grids_pack_size(size_t num_grids, void *grids_, MPI_Comm comm)
struct dist_merge_vtable dist_merge_vtable_grid
static void couple_config_sync_config_file(struct yac_couple_config *couple_config, MPI_Comm comm)
void yac_couple_config_get_field_names(struct yac_couple_config *couple_config, size_t couple_idx, size_t field_couple_idx, char const **src_field_name, const char **tgt_field_name)
struct yac_couple_config * yac_couple_config_new()
static size_t yac_couple_config_get_component_pack_size(struct yac_couple_config_component *component, MPI_Comm comm)
char const * yac_couple_config_get_weight_file_name(struct yac_couple_config *couple_config, size_t couple_idx, size_t field_couple_idx)
struct dist_merge_vtable dist_merge_vtable_couple
static void check_grid_idx(struct yac_couple_config *couple_config, size_t grid_idx, char const *routine_name, int line)
static size_t yac_couple_config_component_add_field_(struct yac_couple_config *couple_config, size_t comp_idx, size_t grid_idx, char const *name, char const *timestep, size_t collection_size)
static void yac_couple_config_pack_component(struct yac_couple_config_component *component, void *buffer, int buffer_size, int *position, MPI_Comm comm)
void yac_couple_config_get_field_couple_component_names(struct yac_couple_config *couple_config, size_t couple_idx, size_t field_couple_idx, char const **src_component_name, char const **tgt_component_name)
static void yac_couple_config_field_free(void *field_)
char const * yac_couple_config_get_component_name(struct yac_couple_config *couple_config, size_t component_idx)
static void yac_couple_config_field_couple_merge(void *a_, void *b_, MPI_Comm comm)
static int yac_couple_config_field_couple_compare(void const *a_, void const *b_)
char const * yac_couple_config_get_source_timestep(struct yac_couple_config *couple_config, size_t couple_idx, size_t field_couple_idx)
static void yac_couple_config_component_free(void *component_)
static void yac_couple_config_pack_string(char const *string, void *buffer, int buffer_size, int *position, MPI_Comm comm)
static void yac_couple_config_pack_grid(struct yac_couple_config_grid *grid, void *buffer, int buffer_size, int *position, MPI_Comm comm)
void yac_couple_config_add_component(struct yac_couple_config *couple_config, char const *name)
size_t yac_couple_config_get_num_grids(struct yac_couple_config *couple_config)
struct dist_merge_vtable dist_merge_vtable_field_couple
char const * yac_couple_config_get_field_grid_name(struct yac_couple_config *couple_config, size_t component_idx, size_t field_idx)
char const * yac_couple_config_get_start_datetime(struct yac_couple_config *couple_config)
static void check_field_idx(struct yac_couple_config *couple_config, size_t component_idx, size_t field_idx, char const *routine_name, int line)
size_t yac_couple_config_get_collection_size(struct yac_couple_config *couple_config, char const *component_name, char const *grid_name, char const *field_name)
static void check_component_idx(struct yac_couple_config *couple_config, size_t component_idx, char const *routine_name, int line)
static void merge_couples(struct yac_couple_config *couple_config, MPI_Comm comm)
static void yac_couple_config_unpack_components(void *buffer, int buffer_size, int *position, size_t *num_components, void *components_, MPI_Comm comm)
char const * yac_couple_config_get_target_timestep(struct yac_couple_config *couple_config, size_t couple_idx, size_t field_couple_idx)
static void yac_couple_config_field_merge(void *a_, void *b_, MPI_Comm comm)
void yac_couple_config_delete(struct yac_couple_config *couple_config)
static void yac_couple_config_pack_field(struct yac_couple_config_field *field, void *buffer, int buffer_size, int *position, MPI_Comm comm)
size_t yac_couple_config_get_num_components(struct yac_couple_config *couple_config)
yac_reduction_type
@ TIME_NONE
@ TIME_ACCUMULATE
@ TIME_MAXIMUM
@ TIME_MINIMUM
@ TIME_AVERAGE
yac_config_filetype
@ YAC_CONFIG_FILETYPE_YAML
YAML format.
@ YAC_CONFIG_FILETYPE_INVALID
invalid format
@ YAC_CONFIG_FILETYPE_JSON
JSON format.
int yac_interp_stack_config_compare(void const *a_, void const *b_)
struct yac_interp_stack_config * yac_interp_stack_config_unpack(void *buffer, int buffer_size, int *position, MPI_Comm comm)
void yac_interp_stack_config_delete(struct yac_interp_stack_config *interp_stack_config)
struct yac_interp_stack_config * yac_interp_stack_config_copy(struct yac_interp_stack_config *interp_stack)
size_t yac_interp_stack_config_get_pack_size(struct yac_interp_stack_config *interp_stack, MPI_Comm comm)
void yac_interp_stack_config_pack(struct yac_interp_stack_config *interp_stack, void *buffer, int buffer_size, int *position, MPI_Comm comm)
double const YAC_FRAC_MASK_UNDEF
double const YAC_FRAC_MASK_NO_VALUE
#define xrealloc(ptr, size)
Definition ppm_xfuncs.h:67
#define xmalloc(size)
Definition ppm_xfuncs.h:66
void yac_qsort_index(void *a_, size_t count, size_t size, int(*compare)(void const *, void const *), size_t *idx)
Definition quicksort.c:185
void(* pack)(size_t, void *, void *, int, int *, MPI_Comm)
size_t(* get_pack_size)(size_t, void *, MPI_Comm)
void(* free_data)(void *)
void(* unpack)(void *, int, int *, size_t *, void *, MPI_Comm)
int(* compare)(void const *, void const *)
void(* merge)(void *, void *, MPI_Comm)
struct yac_couple_config_field * fields
struct yac_couple_config_field_couple * field_couples
struct yac_couple_config_field_couple::@59 target
struct yac_interp_stack_config * interp_stack
struct yac_couple_config_field_couple::@59 source
enum yac_reduction_type coupling_period_operation
struct yac_couple_config_couple * couples
char const * end_datetime
struct yac_couple_config_grid * grids
char const * config_filename
char const * start_datetime
struct yac_couple_config_component * components
enum yac_config_filetype config_filetype
static struct user_input_data_component ** components
Definition yac.c:117
size_t num_grids
Definition yac.c:115
static size_t num_components
Definition yac.c:118
int const YAC_EXCHANGE_TYPE_SOURCE
Definition yac.c:32
int const YAC_EXCHANGE_TYPE_NONE
Definition yac.c:31
int const YAC_EXCHANGE_TYPE_TARGET
Definition yac.c:33
void yac_cdef_calendar(int calendar)
Definition yac.c:616
struct yac_basic_grid ** grids
Definition yac.c:114
#define YAC_MAX_CHARLEN
Definition yac.h:74
#define YAC_ASSERT_F(exp, format,...)
Definition yac_assert.h:19
#define YAC_ASSERT(exp, msg)
Definition yac_assert.h:16
#define yac_mpi_call(call, comm)