YAC 3.21.0
Yet Another Coupler
Loading...
Searching...
No Matches
test_interp_method_avg_parallel.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#include <stdlib.h>
6#include <string.h>
7
8#include "tests.h"
9#include "test_common.h"
10#include "geometry.h"
11#include "interpolation/methods/interp_method.h"
12#include "interpolation/methods/interp_method_avg.h"
13#include "interpolation/methods/interp_method_fixed.h"
14#include "dist_grid_utils.h"
15#include "clipping.h"
16#include "yac_mpi.h"
17
18#include <mpi.h>
19#include <yaxt.h>
20
28size_t num_reorder_types = sizeof(reorder_types) / sizeof(reorder_types[0]);
29static char const * grid_names[2] = {"src_grid", "tgt_grid"};
30
31int main(void) {
32
33 MPI_Init(NULL, NULL);
34
35 xt_initialize(MPI_COMM_WORLD);
36
37 int comm_rank, comm_size;
38 MPI_Comm_rank(MPI_COMM_WORLD, &comm_rank);
39 MPI_Comm_size(MPI_COMM_WORLD, &comm_size);
40 MPI_Barrier(MPI_COMM_WORLD);
41
42 if (comm_size != 3) {
43 PUT_ERR("ERROR: wrong number of processes");
44 xt_finalize();
45 MPI_Finalize();
46 return TEST_EXIT_CODE;
47 }
48
49 MPI_Comm split_comm;
50 MPI_Comm_split(
51 MPI_COMM_WORLD, comm_rank < 2, 0, &split_comm);
52
53 int split_comm_rank, split_comm_size;
54 MPI_Comm_rank(split_comm, &split_comm_rank);
55 MPI_Comm_size(split_comm, &split_comm_size);
56
57 {// Test 1a
58 // linear point interpolation (fixed is backup) with two source processes and a
59 // single target process the source fields have two struct points per grid
60 // (the test checks how the interpolation methods handle NULL target coordinates )
61
62 // the global grid is a 4x2 grid:
63 // 08--14--09--15--10--16--11
64 // | | | |
65 // 08 03 10 04 12 05 13
66 // | | | |
67 // 04--07--05--09--06--11--07
68 // | | | |
69 // 01 00 03 01 05 02 06
70 // | | | |
71 // 00--00--01--02--02--04--03
72 //
73 //---------------
74 // setup
75 //---------------
76
77 int is_tgt = split_comm_size == 1;
78 double coordinates_x[2][4] = {{0.0,1.0,2.0,3.0}, {0.5,1.5,2.5}};
79 double coordinates_y[2][3] = {{0.0,1.0,2.0}, {0.5,1.5,2.5}};
80 size_t const num_cells[2][2] = {{3,2}, {2,2}};
81 size_t local_start[2][2][2] = {{{0,0},{1,0}}, {{0,0}}};
82 size_t local_count[2][2][2] = {{{2,2},{2,2}}, {{2,2}}};
83 int with_halo = 0;
84 for (size_t i = 0; i <= num_cells[is_tgt][0]; ++i)
85 coordinates_x[is_tgt][i] *= YAC_RAD;
86 for (size_t i = 0; i <= num_cells[is_tgt][1]; ++i)
87 coordinates_y[is_tgt][i] *= YAC_RAD;
88
89 struct yac_basic_grid_data grid_data =
90 utest_generate_basic_grid_data_reg2d(
91 coordinates_x[is_tgt], coordinates_y[is_tgt], num_cells[is_tgt],
92 local_start[is_tgt][split_comm_rank],
93 local_count[is_tgt][split_comm_rank], with_halo);
94
95 struct yac_basic_grid * grids[2] =
96 {yac_basic_grid_new(grid_names[is_tgt], grid_data),
97 yac_basic_grid_empty_new(grid_names[is_tgt^1])};
98
99 struct yac_dist_grid_pair * grid_pair =
100 yac_dist_grid_pair_new(grids[0], grids[1], MPI_COMM_WORLD);
101
102 struct yac_interp_field src_fields[] =
103 {{.location = YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX}};
104 size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
105 struct yac_interp_field tgt_field =
106 {.location = YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX};
107
108 struct yac_interp_grid * interp_grid =
109 yac_interp_grid_new(grid_pair, grid_names[0], grid_names[1],
110 num_src_fields, src_fields, tgt_field);
111
112 for (int partial_coverage = 0; partial_coverage <= 1; ++partial_coverage) {
113
114 struct interp_method * method_stack[] =
115 {yac_interp_method_avg_new(
116 YAC_INTERP_AVG_ARITHMETIC, partial_coverage),
117 yac_interp_method_fixed_new(-1.0), NULL};
118
119 struct yac_interp_weights * weights =
120 yac_interp_method_do_search(method_stack, interp_grid);
121
122 for (size_t i = 0; i < num_reorder_types; ++i) {
123
124 struct yac_interpolation * interpolation =
125 yac_interp_weights_get_interpolation(
126 weights, reorder_types[i], 1,
127 YAC_FRAC_MASK_NO_VALUE, 1.0, 0.0, NULL, 1, 1);
128
129 // check generated interpolation
130 {
131 double * src_field = NULL;
132 double ** src_fields = &src_field;
133 double * tgt_field = NULL;
134 double * ref_tgt_field = NULL;
135 double ref_global_tgt_field[9] =
136 {2.5, 3.5, 4.5, 6.5, 7.5, 8.5, -1, -1, -1};
137
138 if (is_tgt) {
139 tgt_field = xmalloc(grid_data.num_vertices * sizeof(*tgt_field));
140 ref_tgt_field = xmalloc(grid_data.num_vertices * sizeof(*ref_tgt_field));
141 for (size_t i = 0; i < grid_data.num_vertices; ++i)
142 ref_tgt_field[i] = ref_global_tgt_field[grid_data.vertex_ids[i]];
143 } else {
144 src_field = xmalloc(grid_data.num_vertices * sizeof(*src_field));
145 for (size_t i = 0; i < grid_data.num_vertices; ++i)
146 src_field[i] = (double)(grid_data.vertex_ids[i]);
147 }
148
149 yac_interpolation_execute(interpolation, &src_fields, &tgt_field);
150
151 if (is_tgt)
152 for (size_t i = 0; i < grid_data.num_vertices; ++i)
153 if (fabs(ref_tgt_field[i] - tgt_field[i]) > 1e-9)
154 PUT_ERR("wrong interpolation result");
155
156 free(src_field);
157 free(tgt_field);
158 free(ref_tgt_field);
159 }
160
161 yac_interpolation_delete(interpolation);
162 }
163
164 yac_interp_weights_delete(weights);
165 yac_interp_method_delete(method_stack);
166 } // partial_coverage
167 yac_interp_grid_delete(interp_grid);
168 yac_dist_grid_pair_delete(grid_pair);
169 yac_basic_grid_delete(grids[1]);
170 yac_basic_grid_delete(grids[0]);
171 }
172
173 {// Test 2a
174 // linear point interpolation (fixed is backup) with two source processes and a
175 // single target process the source fields have two struct points per grid
176 // (the test check how the interpolation methods handle NULL target coordinates )
177
178 // the global grid is a 4x2 grid:
179 // 08--14--09--15--10--16--11
180 // | | | |
181 // 08 03 10 04 12 05 13
182 // | | | |
183 // 04--07--05--09--06--11--07
184 // | | | |
185 // 01 00 03 01 05 02 06
186 // | | | |
187 // 00--00--01--02--02--04--03
188 //
189 // source mask:
190 // 0-------1-------1-------1
191 // | | | |
192 // | | | |
193 // | | | |
194 // 1-------1-------1-------1
195 // | | | |
196 // | | | |
197 // | | | |
198 // 1-------1-------1-------0
199 //
200 //---------------
201 // setup
202 //---------------
203
204 int is_tgt = split_comm_size == 1;
205 double coordinates_x[2][4] = {{0.0,1.0,2.0,3.0}, {0.5,1.5,2.5}};
206 double coordinates_y[2][3] = {{0.0,1.0,2.0}, {0.5,1.5,2.5}};
207 size_t const num_cells[2][2] = {{3,2}, {2,2}};
208 size_t local_start[2][2][2] = {{{0,0},{1,0}}, {{0,0}}};
209 size_t local_count[2][2][2] = {{{2,2},{2,2}}, {{2,2}}};
210 int global_src_mask[12] = {1, 1, 1, 0,
211 1, 1, 1, 1,
212 0, 1, 1, 1};
213 int with_halo = 0;
214 for (size_t i = 0; i <= num_cells[is_tgt][0]; ++i)
215 coordinates_x[is_tgt][i] *= YAC_RAD;
216 for (size_t i = 0; i <= num_cells[is_tgt][1]; ++i)
217 coordinates_y[is_tgt][i] *= YAC_RAD;
218
219 struct yac_basic_grid_data grid_data =
220 utest_generate_basic_grid_data_reg2d(
221 coordinates_x[is_tgt], coordinates_y[is_tgt], num_cells[is_tgt],
222 local_start[is_tgt][split_comm_rank],
223 local_count[is_tgt][split_comm_rank], with_halo);
224
225 struct yac_basic_grid * grids[2] =
226 {yac_basic_grid_new(grid_names[is_tgt], grid_data),
227 yac_basic_grid_empty_new(grid_names[is_tgt^1])};
228 int * src_mask = NULL;
229 if (!is_tgt) {
230 src_mask = xmalloc(grid_data.num_vertices * sizeof(*src_mask));
231 for (size_t i = 0; i < grid_data.num_vertices; ++i)
232 src_mask[i] = global_src_mask[grid_data.vertex_ids[i]];
233 yac_basic_grid_add_mask_nocpy(grids[0], YAC_LOC_CORNER, src_mask, NULL);
234 }
235
236 struct yac_dist_grid_pair * grid_pair =
237 yac_dist_grid_pair_new(grids[0], grids[1], MPI_COMM_WORLD);
238 struct yac_interp_field src_fields[] =
239 {{.location = YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX, .masks_idx = 0}};
240 size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
241 struct yac_interp_field tgt_field =
242 {.location = YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX};
243
244 struct yac_interp_grid * interp_grid =
245 yac_interp_grid_new(grid_pair, grid_names[0], grid_names[1],
246 num_src_fields, src_fields, tgt_field);
247
248 struct interp_method * method_stack[] =
249 {yac_interp_method_avg_new(YAC_INTERP_AVG_ARITHMETIC, 0),
250 yac_interp_method_fixed_new(-1.0), NULL};
251
252 struct yac_interp_weights * weights =
253 yac_interp_method_do_search(method_stack, interp_grid);
254
255 for (size_t i = 0; i < num_reorder_types; ++i) {
256
257 struct yac_interpolation * interpolation =
258 yac_interp_weights_get_interpolation(
259 weights, reorder_types[i], 1,
260 YAC_FRAC_MASK_NO_VALUE, 1.0, 0.0, NULL, 1, 1);
261
262 // check generated interpolation
263 {
264 double * src_field = NULL;
265 double ** src_fields = &src_field;
266 double * tgt_field = NULL;
267 double * ref_tgt_field = NULL;
268 double ref_global_tgt_field[9] =
269 {2.5, 3.5, -1, -1, 7.5, 8.5, -1, -1, -1};
270
271 if (is_tgt) {
272 tgt_field = xmalloc(grid_data.num_vertices * sizeof(*tgt_field));
273 ref_tgt_field = xmalloc(grid_data.num_vertices * sizeof(*ref_tgt_field));
274 for (size_t i = 0; i < grid_data.num_vertices; ++i)
275 ref_tgt_field[i] = ref_global_tgt_field[grid_data.vertex_ids[i]];
276 } else {
277 src_field = xmalloc(grid_data.num_vertices * sizeof(*src_field));
278 for (size_t i = 0; i < grid_data.num_vertices; ++i)
279 src_field[i] = (double)(grid_data.vertex_ids[i]);
280 }
281
282 yac_interpolation_execute(interpolation, &src_fields, &tgt_field);
283
284 if (is_tgt)
285 for (size_t i = 0; i < grid_data.num_vertices; ++i)
286 if (fabs(ref_tgt_field[i] - tgt_field[i]) > 1e-9)
287 PUT_ERR("wrong interpolation result");
288
289 free(src_field);
290 free(tgt_field);
291 free(ref_tgt_field);
292 }
293
294 yac_interpolation_delete(interpolation);
295 }
296
297 yac_interp_weights_delete(weights);
298 yac_interp_method_delete(method_stack);
299 yac_interp_grid_delete(interp_grid);
300 yac_dist_grid_pair_delete(grid_pair);
301 yac_basic_grid_delete(grids[1]);
302 yac_basic_grid_delete(grids[0]);
303 }
304
305 {// Test 3
306 // linear distance weighted (no partial coverage) point interpolation
307 // with two source processes and a single target process
308
309 // the global grid is a 3x3 grid:
310 // 12--21--13--22--14--23--15
311 // | | | |
312 // 15 06 17 07 19 08 20
313 // | | | |
314 // 08--14--09--16--10--18--11
315 // | | | |
316 // 08 03 10 04 12 05 13
317 // | | | |
318 // 04--07--05--09--06--11--07
319 // | | | |
320 // 01 00 03 01 05 02 06
321 // | | | |
322 // 00--00--01--02--02--04--03
323 //
324 // source mask:
325 // 0-------1-------1-------0
326 // | | | |
327 // | | | |
328 // | | | |
329 // 1-------1-------1-------1
330 // | | | |
331 // | | | |
332 // | | | |
333 // 1-------1-------1-------1
334 // | | | |
335 // | | | |
336 // | | | |
337 // 0-------1-------1-------0
338 //
339 //---------------
340 // setup
341 //---------------
342
343 int is_tgt = split_comm_size == 1;
344 double coordinates_x[2][4] = {{0.0,1.0,2.0,3.0}, {0.0,1.0,2.0,3.0}};
345 double coordinates_y[2][4] = {{0.0,1.0,2.0,3.0}, {0.0,1.0,2.0,3.0}};
346 size_t const num_cells[2][2] = {{3,3}, {3,3}};
347 size_t local_start[2][2][2] = {{{0,0},{1,0}}, {{0,0}}};
348 size_t local_count[2][2][2] = {{{2,3},{2,3}}, {{3,3}}};
349 int global_src_mask[16] = {0, 1, 1, 0,
350 1, 1, 1, 1,
351 1, 1, 1, 1,
352 0, 1, 1, 0};
353 int with_halo = 0;
354 for (size_t i = 0; i <= num_cells[is_tgt][0]; ++i)
355 coordinates_x[is_tgt][i] *= YAC_RAD;
356 for (size_t i = 0; i <= num_cells[is_tgt][1]; ++i)
357 coordinates_y[is_tgt][i] *= YAC_RAD;
358
359 struct yac_basic_grid_data grid_data =
360 utest_generate_basic_grid_data_reg2d(
361 coordinates_x[is_tgt], coordinates_y[is_tgt], num_cells[is_tgt],
362 local_start[is_tgt][split_comm_rank],
363 local_count[is_tgt][split_comm_rank], with_halo);
364
365 struct yac_basic_grid * grids[2] =
366 {yac_basic_grid_new(grid_names[is_tgt], grid_data),
367 yac_basic_grid_empty_new(grid_names[is_tgt^1])};
368
369 int * src_mask = NULL;
370 if (!is_tgt) {
371 src_mask = xmalloc(grid_data.num_vertices * sizeof(*src_mask));
372 for (size_t i = 0; i < grid_data.num_vertices; ++i)
373 src_mask[i] = global_src_mask[grid_data.vertex_ids[i]];
374 yac_basic_grid_add_mask_nocpy(grids[0], YAC_LOC_CORNER, src_mask, NULL);
375 }
376
377 yac_coordinate_pointer tgt_cell_coordinates = NULL;
378 if (is_tgt) {
379 double cell_coordinates_x[3] = {0.75,1.75,2.75};
380 double cell_coordinates_y[3] = {0.75,1.75,2.75};
381 tgt_cell_coordinates = xmalloc(9 * sizeof(*tgt_cell_coordinates));
382 for (int i = 0, k = 0; i < 3; ++i)
383 for (int j = 0; j < 3; ++j, ++k)
384 LLtoXYZ_deg(
385 cell_coordinates_x[j], cell_coordinates_y[i],
386 tgt_cell_coordinates[k]);
387 yac_basic_grid_add_coordinates_nocpy(grids[0], YAC_LOC_CELL, tgt_cell_coordinates);
388 }
389
390 struct yac_dist_grid_pair * grid_pair =
391 yac_dist_grid_pair_new(grids[0], grids[1], MPI_COMM_WORLD);
392
393 struct yac_interp_field src_fields[] =
394 {{.location = YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX, .masks_idx = 0}};
395 size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
396 struct yac_interp_field tgt_field =
397 {.location = YAC_LOC_CELL, .coordinates_idx = 0, .masks_idx = SIZE_MAX};
398
399 struct yac_interp_grid * interp_grid =
400 yac_interp_grid_new(grid_pair, grid_names[0], grid_names[1],
401 num_src_fields, src_fields, tgt_field);
402
403 struct interp_method * method_stack[] =
404 {yac_interp_method_avg_new(YAC_INTERP_AVG_DIST, 0),
405 yac_interp_method_fixed_new(-1.0), NULL};
406
407 struct yac_interp_weights * weights =
408 yac_interp_method_do_search(method_stack, interp_grid);
409
410 for (size_t i = 0; i < num_reorder_types; ++i) {
411
412 struct yac_interpolation * interpolation =
413 yac_interp_weights_get_interpolation(
414 weights, reorder_types[i], 1,
415 YAC_FRAC_MASK_NO_VALUE, 1.0, 0.0, NULL, 1, 1);
416
417 // check generated interpolation
418 {
419 double * src_field = NULL;
420 double ** src_fields = &src_field;
421 double * tgt_field = NULL;
422 double * ref_tgt_field = NULL;
423 double ref_global_tgt_field[9] =
424 {-1.0,//0.14963* 0.0 + 0.20075*( 1.0+ 4.0) + 0.44888* 5.0,
425 0.14963* 1.0 + 0.20075*( 2.0+ 5.0) + 0.44888* 6.0,
426 -1.0,//0.14963* 2.0 + 0.20075*( 3.0+ 6.0) + 0.44888* 7.0,
427 0.14963* 4.0 + 0.20075*( 5.0+ 8.0) + 0.44888* 9.0,
428 0.14963* 5.0 + 0.20075*( 6.0+ 9.0) + 0.44888*10.0,
429 0.14963* 6.0 + 0.20075*( 7.0+10.0) + 0.44888*11.0,
430 -1.0,//0.14963* 8.0 + 0.20075*( 9.0+12.0) + 0.44888*13.0,
431 0.14963* 9.0 + 0.20075*(10.0+13.0) + 0.44888*14.0,
432 -1.0};//0.14963*10.0 + 0.20075*(11.0+14.0) + 0.44888*15.0};
433
434 if (is_tgt) {
435 tgt_field = xmalloc(grid_data.num_vertices * sizeof(*tgt_field));
436 ref_tgt_field = xmalloc(grid_data.num_vertices * sizeof(*ref_tgt_field));
437 for (size_t i = 0; i < grid_data.num_vertices; ++i)
438 ref_tgt_field[i] = ref_global_tgt_field[grid_data.vertex_ids[i]];
439 } else {
440 src_field = xmalloc(grid_data.num_vertices * sizeof(*src_field));
441 for (size_t i = 0; i < grid_data.num_vertices; ++i)
442 src_field[i] = (double)(grid_data.vertex_ids[i]);
443 }
444
445 yac_interpolation_execute(interpolation, &src_fields, &tgt_field);
446
447 if (is_tgt)
448 for (size_t i = 0; i < grid_data.num_cells; ++i)
449 if (fabs(ref_tgt_field[i] - tgt_field[i]) > 1e-3)
450 PUT_ERR("wrong interpolation result");
451
452 free(src_field);
453 free(tgt_field);
454 free(ref_tgt_field);
455 }
456
457 yac_interpolation_delete(interpolation);
458 }
459
460 yac_interp_weights_delete(weights);
461 yac_interp_method_delete(method_stack);
462 yac_interp_grid_delete(interp_grid);
463 yac_dist_grid_pair_delete(grid_pair);
464 yac_basic_grid_delete(grids[1]);
465 yac_basic_grid_delete(grids[0]);
466 }
467
468 {// Test 4
469 // linear distance weighted (with partial coverage) point interpolation
470 // with two source processes and a single target process
471
472 // the global grid is a 3x3 grid:
473 // 12--21--13--22--14--23--15
474 // | | | |
475 // 15 06 17 07 19 08 20
476 // | | | |
477 // 08--14--09--16--10--18--11
478 // | | | |
479 // 08 03 10 04 12 05 13
480 // | | | |
481 // 04--07--05--09--06--11--07
482 // | | | |
483 // 01 00 03 01 05 02 06
484 // | | | |
485 // 00--00--01--02--02--04--03
486 //
487 // source mask:
488 // 0-------1-------1-------0
489 // | | | |
490 // | | | |
491 // | | | |
492 // 1-------0-------1-------1
493 // | | | |
494 // | | | |
495 // | | | |
496 // 1-------1-------0-------0
497 // | | | |
498 // | | | |
499 // | | | |
500 // 0-------1-------0-------0
501 //
502 //---------------
503 // setup
504 //---------------
505
506 int is_tgt = split_comm_size == 1;
507 double coordinates_x[2][4] = {{0.0,1.0,2.0,3.0}, {0.0,1.0,2.0,3.0}};
508 double coordinates_y[2][4] = {{0.0,1.0,2.0,3.0}, {0.0,1.0,2.0,3.0}};
509 size_t const num_cells[2][2] = {{3,3}, {3,3}};
510 size_t local_start[2][2][2] = {{{0,0},{1,0}}, {{0,0}}};
511 size_t local_count[2][2][2] = {{{2,3},{2,3}}, {{3,3}}};
512 int global_src_mask[16] = {0, 1, 0, 0,
513 1, 1, 0, 0,
514 1, 0, 1, 1,
515 0, 1, 1, 0};
516 int with_halo = 0;
517 for (size_t i = 0; i <= num_cells[is_tgt][0]; ++i)
518 coordinates_x[is_tgt][i] *= YAC_RAD;
519 for (size_t i = 0; i <= num_cells[is_tgt][1]; ++i)
520 coordinates_y[is_tgt][i] *= YAC_RAD;
521
522 struct yac_basic_grid_data grid_data =
523 utest_generate_basic_grid_data_reg2d(
524 coordinates_x[is_tgt], coordinates_y[is_tgt], num_cells[is_tgt],
525 local_start[is_tgt][split_comm_rank],
526 local_count[is_tgt][split_comm_rank], with_halo);
527
528 struct yac_basic_grid * grids[2] =
529 {yac_basic_grid_new(grid_names[is_tgt], grid_data),
530 yac_basic_grid_empty_new(grid_names[is_tgt^1])};
531
532 int * src_mask = NULL;
533 if (!is_tgt) {
534 src_mask = xmalloc(grid_data.num_vertices * sizeof(*src_mask));
535 for (size_t i = 0; i < grid_data.num_vertices; ++i)
536 src_mask[i] = global_src_mask[grid_data.vertex_ids[i]];
537 yac_basic_grid_add_mask_nocpy(grids[0], YAC_LOC_CORNER, src_mask, NULL);
538 }
539
540 yac_coordinate_pointer tgt_cell_coordinates = NULL;
541 if (is_tgt) {
542 double cell_coordinates_x[3] = {0.75,1.75,2.75};
543 double cell_coordinates_y[3] = {0.75,1.75,2.75};
544 tgt_cell_coordinates = xmalloc(9 * sizeof(*tgt_cell_coordinates));
545 for (int i = 0, k = 0; i < 3; ++i)
546 for (int j = 0; j < 3; ++j, ++k)
547 LLtoXYZ_deg(
548 cell_coordinates_x[j], cell_coordinates_y[i],
549 tgt_cell_coordinates[k]);
550 yac_basic_grid_add_coordinates_nocpy(grids[0], YAC_LOC_CELL, tgt_cell_coordinates);
551 }
552
553 struct yac_dist_grid_pair * grid_pair =
554 yac_dist_grid_pair_new(grids[0], grids[1], MPI_COMM_WORLD);
555
556 struct yac_interp_field src_fields[] =
557 {{.location = YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX, .masks_idx = 0}};
558 size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
559 struct yac_interp_field tgt_field =
560 {.location = YAC_LOC_CELL, .coordinates_idx = 0, .masks_idx = SIZE_MAX};
561
562 struct yac_interp_grid * interp_grid =
563 yac_interp_grid_new(grid_pair, grid_names[0], grid_names[1],
564 num_src_fields, src_fields, tgt_field);
565
566 struct interp_method * method_stack[] =
567 {yac_interp_method_avg_new(YAC_INTERP_AVG_DIST, 1),
568 yac_interp_method_fixed_new(-1.0), NULL};
569
570 struct yac_interp_weights * weights =
571 yac_interp_method_do_search(method_stack, interp_grid);
572
573 for (size_t i = 0; i < num_reorder_types; ++i) {
574
575 struct yac_interpolation * interpolation =
576 yac_interp_weights_get_interpolation(
577 weights, reorder_types[i], 1,
578 YAC_FRAC_MASK_NO_VALUE, 1.0, 0.0, NULL, 1, 1);
579
580 // check generated interpolation
581 {
582 double * src_field = NULL;
583 double ** src_fields = &src_field;
584 double * tgt_field = NULL;
585 double * ref_tgt_field = NULL;
586 double inv_dist_a = sqrt(8.0);
587 double inv_dist_b = sqrt(8.0/5.0);
588 double inv_dist_c = sqrt(8.0/9.0);
589 double ref_global_tgt_field[9] =
590 {(inv_dist_b*( 1.0+ 4.0) + inv_dist_a* 5.0)/(inv_dist_a + 2.0 * inv_dist_b),
591 (inv_dist_c* 1.0 + inv_dist_b* 5.0)/(inv_dist_b + inv_dist_c),
592 -1.0,
593 (inv_dist_c* 4.0 + inv_dist_b*( 5.0+ 8.0))/(2.0 * inv_dist_b + inv_dist_c),
594 (inv_dist_c* 5.0 + inv_dist_a*10.0)/(inv_dist_a + inv_dist_c),
595 (inv_dist_b*10.0 + inv_dist_a*11.0)/(inv_dist_a + inv_dist_b),
596 (inv_dist_c* 8.0 + inv_dist_a*13.0)/(inv_dist_a + inv_dist_c),
597 (inv_dist_b*(10.0+13.0) + inv_dist_a*14.0)/(inv_dist_a + 2.0 * inv_dist_b),
598 (inv_dist_c*10.0 + inv_dist_b*(11.0+14.0))/(2.0 * inv_dist_b + inv_dist_c)};
599
600 if (is_tgt) {
601 tgt_field = xmalloc(grid_data.num_vertices * sizeof(*tgt_field));
602 ref_tgt_field = xmalloc(grid_data.num_vertices * sizeof(*ref_tgt_field));
603 for (size_t i = 0; i < grid_data.num_vertices; ++i)
604 ref_tgt_field[i] = ref_global_tgt_field[grid_data.vertex_ids[i]];
605 } else {
606 src_field = xmalloc(grid_data.num_vertices * sizeof(*src_field));
607 for (size_t i = 0; i < grid_data.num_vertices; ++i)
608 src_field[i] = (double)(grid_data.vertex_ids[i]);
609 }
610
611 yac_interpolation_execute(interpolation, &src_fields, &tgt_field);
612
613 if (is_tgt)
614 for (size_t i = 0; i < grid_data.num_cells; ++i)
615 if (fabs(ref_tgt_field[i] - tgt_field[i]) > 1e-3)
616 PUT_ERR("wrong interpolation result");
617
618 free(src_field);
619 free(tgt_field);
620 free(ref_tgt_field);
621 }
622
623 yac_interpolation_delete(interpolation);
624 }
625
626 yac_interp_weights_delete(weights);
627 yac_interp_method_delete(method_stack);
628 yac_interp_grid_delete(interp_grid);
629 yac_dist_grid_pair_delete(grid_pair);
630 yac_basic_grid_delete(grids[1]);
631 yac_basic_grid_delete(grids[0]);
632 }
633
634 {// Test 5
635 // linear distance weighted (with partial coverage) point interpolation
636 // with two source processes and a single target process
637 // (target points are on vertices of the source grid
638
639 // the global grid is a 3x3 grid:
640 // 12--21--13--22--14--23--15
641 // | | | |
642 // 15 06 17 07 19 08 20
643 // | | | |
644 // 08--14--09--16--10--18--11
645 // | | | |
646 // 08 03 10 04 12 05 13
647 // | | | |
648 // 04--07--05--09--06--11--07
649 // | | | |
650 // 01 00 03 01 05 02 06
651 // | | | |
652 // 00--00--01--02--02--04--03
653 //
654 //---------------
655 // setup
656 //---------------
657
658 int is_tgt = split_comm_size == 1;
659 double coordinates_x[2][4] = {{0.0,1.0,2.0,3.0}, {0.0,1.0,2.0,3.0}};
660 double coordinates_y[2][4] = {{0.0,1.0,2.0,3.0}, {0.0,1.0,2.0,3.0}};
661 size_t const num_cells[2][2] = {{3,3}, {3,3}};
662 size_t local_start[2][2][2] = {{{0,0},{1,0}}, {{0,0}}};
663 size_t local_count[2][2][2] = {{{2,3},{2,3}}, {{3,3}}};
664 int with_halo = 0;
665 for (size_t i = 0; i <= num_cells[is_tgt][0]; ++i)
666 coordinates_x[is_tgt][i] *= YAC_RAD;
667 for (size_t i = 0; i <= num_cells[is_tgt][1]; ++i)
668 coordinates_y[is_tgt][i] *= YAC_RAD;
669
670 struct yac_basic_grid_data grid_data =
671 utest_generate_basic_grid_data_reg2d(
672 coordinates_x[is_tgt], coordinates_y[is_tgt], num_cells[is_tgt],
673 local_start[is_tgt][split_comm_rank],
674 local_count[is_tgt][split_comm_rank], with_halo);
675
676 struct yac_basic_grid * grids[2] =
677 {yac_basic_grid_new(grid_names[is_tgt], grid_data),
678 yac_basic_grid_empty_new(grid_names[is_tgt^1])};
679
680 yac_coordinate_pointer tgt_cell_coordinates = NULL;
681 if (is_tgt) {
682 double cell_coordinates_x[3] = {1.0,2.0,3.0};
683 double cell_coordinates_y[3] = {1.0,2.0,3.0};
684 tgt_cell_coordinates = xmalloc(9 * sizeof(*tgt_cell_coordinates));
685 for (int i = 0, k = 0; i < 3; ++i)
686 for (int j = 0; j < 3; ++j, ++k)
687 LLtoXYZ_deg(
688 cell_coordinates_x[j], cell_coordinates_y[i],
689 tgt_cell_coordinates[k]);
690 yac_basic_grid_add_coordinates_nocpy(grids[0], YAC_LOC_CELL, tgt_cell_coordinates);
691 }
692
693 struct yac_dist_grid_pair * grid_pair =
694 yac_dist_grid_pair_new(grids[0], grids[1], MPI_COMM_WORLD);
695
696 struct yac_interp_field src_fields[] =
697 {{.location = YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX}};
698 size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
699 struct yac_interp_field tgt_field =
700 {.location = YAC_LOC_CELL, .coordinates_idx = 0, .masks_idx = SIZE_MAX};
701
702 struct yac_interp_grid * interp_grid =
703 yac_interp_grid_new(grid_pair, grid_names[0], grid_names[1],
704 num_src_fields, src_fields, tgt_field);
705
706 struct interp_method * method_stack[] =
707 {yac_interp_method_avg_new(YAC_INTERP_AVG_DIST, 1),
708 yac_interp_method_fixed_new(-1.0), NULL};
709
710 struct yac_interp_weights * weights =
711 yac_interp_method_do_search(method_stack, interp_grid);
712
713 for (size_t i = 0; i < num_reorder_types; ++i) {
714
715 struct yac_interpolation * interpolation =
716 yac_interp_weights_get_interpolation(
717 weights, reorder_types[i], 1,
718 YAC_FRAC_MASK_NO_VALUE, 1.0, 0.0, NULL, 1, 1);
719
720 // check generated interpolation
721 {
722 double * src_field = NULL;
723 double ** src_fields = &src_field;
724 double * tgt_field = NULL;
725 double * ref_tgt_field = NULL;
726 double ref_global_tgt_field[9] =
727 {5.0, 6.0, 7.0, 9.0, 10.0, 11.0, 13.0, 14.0, 15.0};
728
729 if (is_tgt) {
730 tgt_field = xmalloc(grid_data.num_vertices * sizeof(*tgt_field));
731 ref_tgt_field = xmalloc(grid_data.num_vertices * sizeof(*ref_tgt_field));
732 for (size_t i = 0; i < grid_data.num_vertices; ++i)
733 ref_tgt_field[i] = ref_global_tgt_field[grid_data.vertex_ids[i]];
734 } else {
735 src_field = xmalloc(grid_data.num_vertices * sizeof(*src_field));
736 for (size_t i = 0; i < grid_data.num_vertices; ++i)
737 src_field[i] = (double)(grid_data.vertex_ids[i]);
738 }
739
740 yac_interpolation_execute(interpolation, &src_fields, &tgt_field);
741
742 if (is_tgt)
743 for (size_t i = 0; i < grid_data.num_cells; ++i)
744 if (fabs(ref_tgt_field[i] - tgt_field[i]) > 1e-9)
745 PUT_ERR("wrong interpolation result");
746
747 free(src_field);
748 free(tgt_field);
749 free(ref_tgt_field);
750 }
751
752 yac_interpolation_delete(interpolation);
753 }
754
755 yac_interp_weights_delete(weights);
756 yac_interp_method_delete(method_stack);
757 yac_interp_grid_delete(interp_grid);
758 yac_dist_grid_pair_delete(grid_pair);
759 yac_basic_grid_delete(grids[1]);
760 yac_basic_grid_delete(grids[0]);
761 }
762
763 {// Test 6
764 // tests cell based source field
765
766 // the global grid is a 7x7 grid:
767 // 56-----57-----58-----59-----60-----61-----62-----63
768 // | | | | | | | |
769 // | 42 | 43 | 44 | 45 | 46 | 47 | 48 |
770 // | | | | | | | |
771 // 48-----49-----50-----51-----52-----53-----54-----55
772 // | | | | | | | |
773 // | 35 | 36 | 37 | 38 | 39 | 40 | 41 |
774 // | | | | | | | |
775 // 40-----41-----42-----43-----44-----45-----46-----47
776 // | | | | | | | |
777 // | 28 | 29 | 30 | 31 | 32 | 33 | 34 |
778 // | | | | | | | |
779 // 32-----33-----34-----35-----36-----37-----38-----39
780 // | | | | | | | |
781 // | 21 | 22 | 23 | 24 | 25 | 26 | 27 |
782 // | | | | | | | |
783 // 24-----25-----26-----27-----28-----29-----30-----31
784 // | | | | | | | |
785 // | 14 | 15 | 16 | 17 | 18 | 19 | 20 |
786 // | | | | | | | |
787 // 16-----17-----18-----19-----20-----21-----22-----23
788 // | | | | | | | |
789 // | 07 | 08 | 09 | 10 | 11 | 12 | 13 |
790 // | | | | | | | |
791 // 08-----09-----10-----11-----12-----13-----14-----15
792 // | | | | | | | |
793 // | 00 | 01 | 02 | 03 | 04 | 05 | 06 |
794 // | | | | | | | |
795 // 00-----01-----02-----03-----04-----05-----06-----07
796 //
797 //---------------
798 // setup
799 //---------------
800
801 int is_tgt = split_comm_size == 1;
802 double coordinates_x[8] = {0.0,1.0,2.0,3.0,4.0,5.0,6.0,7.0};
803 double coordinates_y[8] = {0.0,1.0,2.0,3.0,4.0,5.0,6.0,7.0};
804 size_t const num_cells[2] = {7,7};
805 size_t local_start[2][2][2] = {{{0,0},{0,4}}, {{0,0}}};
806 size_t local_count[2][2][2] = {{{7,4},{7,3}}, {{7,7}}};
807 int with_halo = 0;
808 for (size_t i = 0; i <= num_cells[0]; ++i)
809 coordinates_x[i] *= YAC_RAD;
810 for (size_t i = 0; i <= num_cells[1]; ++i)
811 coordinates_y[i] *= YAC_RAD;
812
813 struct yac_basic_grid_data grid_data =
814 utest_generate_basic_grid_data_reg2d(
815 coordinates_x, coordinates_y, num_cells,
816 local_start[is_tgt][split_comm_rank],
817 local_count[is_tgt][split_comm_rank], with_halo);
818
819 struct yac_basic_grid * grids[2] =
820 {yac_basic_grid_new(grid_names[is_tgt], grid_data),
821 yac_basic_grid_empty_new(grid_names[is_tgt^1])};
822
823 yac_coordinate_pointer src_point_coordinates = NULL;
824 if (!is_tgt) {
825 src_point_coordinates =
826 xmalloc(grid_data.num_cells * sizeof(*src_point_coordinates));
827 for (size_t i = 0; i < grid_data.num_cells; ++i) {
828 double * middle_point = src_point_coordinates[i];
829 for (size_t k = 0; k < 3; ++k) middle_point[k] = 0.0;
830 size_t * curr_vertices =
831 grid_data.cell_to_vertex + grid_data.cell_to_vertex_offsets[i];
832 size_t curr_num_vertices = grid_data.num_vertices_per_cell[i];
833 for (size_t j = 0; j < curr_num_vertices; ++j) {
834 double * curr_vertex_coord =
835 grid_data.vertex_coordinates[curr_vertices[j]];
836 for (size_t k = 0; k < 3; ++k)
837 middle_point[k] += curr_vertex_coord[k];
838 }
839 normalise_vector(middle_point);
840 }
841 yac_basic_grid_add_coordinates_nocpy(grids[0], YAC_LOC_CELL, src_point_coordinates);
842 }
843
844 struct yac_dist_grid_pair * grid_pair =
845 yac_dist_grid_pair_new(grids[0], grids[1], MPI_COMM_WORLD);
846 struct yac_interp_field src_fields[] =
847 {{.location = YAC_LOC_CELL, .coordinates_idx = 0, .masks_idx = SIZE_MAX}};
848 size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
849 struct yac_interp_field tgt_field =
850 {.location = YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX};
851
852 struct yac_interp_grid * interp_grid =
853 yac_interp_grid_new(grid_pair, grid_names[0], grid_names[1],
854 num_src_fields, src_fields, tgt_field);
855
856 struct interp_method * method_stack[] =
857 {yac_interp_method_avg_new(YAC_INTERP_AVG_ARITHMETIC, 0),
858 yac_interp_method_fixed_new(-1.0), NULL};
859
860 struct yac_interp_weights * weights =
861 yac_interp_method_do_search(method_stack, interp_grid);
862
863 for (size_t i = 0; i < num_reorder_types; ++i) {
864
865 struct yac_interpolation * interpolation =
866 yac_interp_weights_get_interpolation(
867 weights, reorder_types[i], 1,
868 YAC_FRAC_MASK_NO_VALUE, 1.0, 0.0, NULL, 1, 1);
869
870 // check generated interpolation
871 {
872 double * src_field = NULL;
873 double ** src_fields = &src_field;
874 double * tgt_field = NULL;
875 double * ref_tgt_field = NULL;
876 double ref_global_tgt_field[64] =
877 {-4,-4,-4,-4,-4,-4,-4,-4,
878 -4, 0+ 1+ 7+ 8, 1+ 2+ 8+ 9, 2+ 3+ 9+10, 3+ 4+10+11, 4+ 5+11+12, 5+ 6+12+13,-4,
879 -4, 7+ 8+14+15, 8+ 9+15+16, 9+10+16+17,10+11+17+18,11+12+18+19,12+13+19+20,-4,
880 -4,14+15+21+22,15+16+22+23,16+17+23+24,17+18+24+25,18+19+25+26,19+20+26+27,-4,
881 -4,21+22+28+29,22+23+29+30,23+24+30+31,24+25+31+32,25+26+32+33,26+27+33+34,-4,
882 -4,28+29+35+36,29+30+36+37,30+31+37+38,31+32+38+39,32+33+39+40,33+34+40+41,-4,
883 -4,35+36+42+43,36+37+43+44,37+38+44+45,38+39+45+46,39+40+46+47,40+41+47+48,-4,
884 -4,-4,-4,-4,-4,-4,-4,-4};
885 for (size_t i = 0;
886 i < sizeof(ref_global_tgt_field) / sizeof(ref_global_tgt_field[0]);
887 ++i)
888 ref_global_tgt_field[i] /= 4.0;
889
890 if (is_tgt) {
891 tgt_field = xmalloc(grid_data.num_vertices * sizeof(*tgt_field));
892 ref_tgt_field = xmalloc(grid_data.num_vertices * sizeof(*ref_tgt_field));
893 for (size_t i = 0; i < grid_data.num_vertices; ++i)
894 ref_tgt_field[i] = ref_global_tgt_field[grid_data.vertex_ids[i]];
895 } else {
896 src_field = xmalloc(grid_data.num_cells * sizeof(*src_field));
897 for (size_t i = 0; i < grid_data.num_cells; ++i)
898 src_field[i] = (double)(grid_data.cell_ids[i]);
899 }
900
901 yac_interpolation_execute(interpolation, &src_fields, &tgt_field);
902
903 if (is_tgt)
904 for (size_t i = 0; i < grid_data.num_vertices; ++i)
905 if (fabs(ref_tgt_field[i] - tgt_field[i]) > 1e-9)
906 PUT_ERR("wrong interpolation result");
907
908 free(src_field);
909 free(tgt_field);
910 free(ref_tgt_field);
911 }
912
913 yac_interpolation_delete(interpolation);
914 }
915
916 yac_interp_weights_delete(weights);
917 yac_interp_method_delete(method_stack);
918 yac_interp_grid_delete(interp_grid);
919 yac_dist_grid_pair_delete(grid_pair);
920 yac_basic_grid_delete(grids[1]);
921 yac_basic_grid_delete(grids[0]);
922 }
923
924 {// Test 7
925 // barycentric coordinates weighted (with partial coverage) point
926 // interpolation with two source processes and a single target process
927
928 // the global source grid is a 1x1 grid:
929 // 02--03--03
930 // | |
931 // 01 00 02
932 // | |
933 // 00--00--01
934 //
935 //---------------
936 // setup
937 //---------------
938
939 int is_tgt = split_comm_size == 1;
940 double src_coordinates_x[2][2] = {{-0.5,0.5}, {0.5,-0.5}};
941 double src_coordinates_y[2][2][2] =
942 {{{-0.5,0.5}, {0.5,-0.5}}, {{89,90}, {90,89}}};
943 double src_coordinates_x_unstruct[4];
944 double src_coordinates_y_unstruct[4];
945 size_t const src_num_cells[2] = {1,1};
946 size_t src_local_start[2] = {0,0};
947 size_t src_local_count[2] = {1,1};
948 int src_cell_mask[4] = {1,1, 0,1};
949 yac_int src_reorder[2][2][4] =
950 {{{0,1,2,3},{2,3,0,1}},{{1,0,3,2},{3,2,1,0}}};
951 double tgt_coordinates_x[5] = {-0.5,-0.25,0.0,0.25,0.5};
952 double tgt_coordinates_y[2][5] =
953 {{-0.5,-0.25,0.0,0.25,0.5}, {89.00,89.25,89.50,89.75,90.00}};
954 size_t const tgt_num_cells[2] = {4,4};
955 size_t tgt_local_start[2] = {0,0};
956 size_t tgt_local_count[2] = {4,4};
957
958 for (size_t i = 0; i < 2 * 2; ++i)
959 (&(src_coordinates_x[0][0]))[i] *= YAC_RAD;
960 for (size_t i = 0; i < 2 * 2 * 2; ++i)
961 (&(src_coordinates_y[0][0][0]))[i] *= YAC_RAD;
962 for (size_t i = 0; i < 5; ++i) tgt_coordinates_x[i] *= YAC_RAD;
963 for (size_t i = 0; i < 2 * 5; ++i)
964 (&(tgt_coordinates_y[0][0]))[i] *= YAC_RAD;
965
966 for (int partial_coverage = 0; partial_coverage < 2; ++partial_coverage) {
967 for (int is_reg2d_grid = 0; is_reg2d_grid < 2; ++is_reg2d_grid) {
968 for (int at_pole = 0; at_pole < 2; ++at_pole) {
969 for (int with_mask = 0; with_mask < 2; ++with_mask) {
970 for (int order_x = 0; order_x < 2; ++order_x) {
971 for (int order_y = 0; order_y < 2; ++order_y) {
972
973 struct yac_basic_grid_data grid_data;
974
975 if (is_tgt) {
976 grid_data =
977 utest_generate_basic_grid_data_reg2d(
978 tgt_coordinates_x, tgt_coordinates_y[at_pole],
979 tgt_num_cells, tgt_local_start, tgt_local_count, 0);
980 } else {
981 if (is_reg2d_grid) {
982 grid_data =
983 utest_generate_basic_grid_data_reg2d(
984 src_coordinates_x[order_x],
985 src_coordinates_y[at_pole][order_y], src_num_cells,
986 src_local_start, src_local_count, 0);
987 } else {
988 for (int i = 0; i < 4; ++i) {
989 src_coordinates_x_unstruct[i] =
990 src_coordinates_x[order_x][i&1];
991 src_coordinates_y_unstruct[i] =
992 src_coordinates_y[at_pole][order_y][i>>1];
993 }
994 grid_data =
995 yac_generate_basic_grid_data_unstruct(
996 4, 1, (int[]){4},
997 src_coordinates_x_unstruct,
998 src_coordinates_y_unstruct, (int[]){0,1,3,2});
999 }
1000 }
1001
1002 if (!is_tgt) {
1003
1004 yac_int * vertex_ids =
1005 xmalloc(grid_data.num_vertices * sizeof(*vertex_ids));
1006 for (size_t i = 0; i < grid_data.num_vertices; ++i)
1007 vertex_ids[i] = (yac_int)(src_reorder[order_x][order_y][i]);
1008 grid_data.vertex_ids = vertex_ids;
1009 }
1010
1011 struct yac_basic_grid * grids[2] =
1012 {yac_basic_grid_new(grid_names[is_tgt], grid_data),
1013 yac_basic_grid_empty_new(grid_names[is_tgt^1])};
1014
1015 if (with_mask && !is_tgt) {
1016 int * mask = xmalloc(grid_data.num_vertices * sizeof(*mask));
1017 for (size_t i = 0; i < grid_data.num_vertices; ++i)
1018 mask[i] =
1019 src_cell_mask[src_reorder[order_x][order_y][i]];
1020 yac_basic_grid_add_mask_nocpy(grids[0], YAC_LOC_CORNER, mask, NULL);
1021 }
1022
1023 struct yac_dist_grid_pair * grid_pair =
1024 yac_dist_grid_pair_new(grids[0], grids[1], MPI_COMM_WORLD);
1025
1026 struct yac_interp_field src_fields[] =
1027 {{.location = YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX,
1028 .masks_idx = (with_mask)?0:SIZE_MAX}};
1029 size_t num_src_fields =
1030 sizeof(src_fields) / sizeof(src_fields[0]);
1031 struct yac_interp_field tgt_field =
1032 {.location = YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX,
1033 .masks_idx = SIZE_MAX};
1034
1035 struct yac_interp_grid * interp_grid =
1036 yac_interp_grid_new(grid_pair, grid_names[0], grid_names[1],
1037 num_src_fields, src_fields, tgt_field);
1038
1039 struct interp_method * method_stack[] =
1040 {yac_interp_method_avg_new(
1041 YAC_INTERP_AVG_BARY, partial_coverage),
1042 yac_interp_method_fixed_new(-1.0), NULL};
1043
1044 struct yac_interp_weights * weights =
1045 yac_interp_method_do_search(method_stack, interp_grid);
1046
1047 for (size_t i = 0; i < num_reorder_types; ++i) {
1048
1049 struct yac_interpolation * interpolation =
1050 yac_interp_weights_get_interpolation(
1051 weights, reorder_types[i], 1,
1052 YAC_FRAC_MASK_NO_VALUE, 1.0, 0.0, NULL, 1, 1);
1053
1054 // check generated interpolation
1055 {
1056 double * src_field = NULL;
1057 double ** src_fields = &src_field;
1058 double * tgt_field = NULL;
1059 double * ref_tgt_field = NULL;
1060 double ref_global_tgt_field[25];
1061 double ref_weights[2][2][2][25][4] =
1062 {{{{{1.00, 0.00, 0.00, 0.00},
1063 {0.75, 0.25, 0.00, 0.00},
1064 {0.50, 0.50, 0.00, 0.00},
1065 {0.25, 0.75, 0.00, 0.00},
1066 {0.00, 1.00, 0.00, 0.00},
1067 {0.75, 0.00, 0.25, 0.00},
1068 {0.75, 0.00, 0.00, 0.25},
1069 {0.50, 0.25, 0.00, 0.25},
1070 {0.25, 0.50, 0.00, 0.25},
1071 {0.00, 0.75, 0.00, 0.25},
1072 {0.50, 0.00, 0.50, 0.00},
1073 {0.50, 0.00, 0.25, 0.25},
1074 {0.50, 0.00, 0.00, 0.50},
1075 {0.25, 0.25, 0.00, 0.50},
1076 {0.00, 0.50, 0.00, 0.50},
1077 {0.25, 0.00, 0.75, 0.00},
1078 {0.25, 0.00, 0.50, 0.25},
1079 {0.25, 0.00, 0.25, 0.50},
1080 {0.25, 0.00, 0.00, 0.75},
1081 {0.00, 0.25, 0.00, 0.75},
1082 {0.00, 0.00, 1.00, 0.00},
1083 {0.00, 0.00, 0.75, 0.25},
1084 {0.00, 0.00, 0.50, 0.50},
1085 {0.00, 0.00, 0.25, 0.75},
1086 {0.00, 0.00, 0.00, 1.00}},
1087 {{1.00, 0.00, 0.00, 0.00},
1088 {0.75, 0.25, 0.00, 0.00},
1089 {0.50, 0.50, 0.00, 0.00},
1090 {0.25, 0.75, 0.00, 0.00},
1091 {0.00, 1.00, 0.00, 0.00},
1092 {0.00, 0.00, 0.00, 0.00},//{0.75, 0.00, 0.25, 0.00}
1093 {0.75, 0.00, 0.00, 0.25},
1094 {0.50, 0.25, 0.00, 0.25},
1095 {0.25, 0.50, 0.00, 0.25},
1096 {0.00, 0.75, 0.00, 0.25},
1097 {0.00, 0.00, 0.00, 0.00},//{0.50, 0.00, 0.50, 0.00}
1098 {0.00, 0.00, 0.00, 0.00},//{0.50, 0.00, 0.25, 0.25}
1099 {0.50, 0.00, 0.00, 0.50},
1100 {0.25, 0.25, 0.00, 0.50},
1101 {0.00, 0.50, 0.00, 0.50},
1102 {0.00, 0.00, 0.00, 0.00},//{0.25, 0.00, 0.75, 0.00}
1103 {0.00, 0.00, 0.00, 0.00},//{0.25, 0.00, 0.50, 0.25}
1104 {0.00, 0.00, 0.00, 0.00},//{0.25, 0.00, 0.25, 0.50}
1105 {0.25, 0.00, 0.00, 0.75},
1106 {0.00, 0.25, 0.00, 0.75},
1107 {0.00, 0.00, 0.00, 0.00},//{0.00, 0.00, 1.00, 0.00}
1108 {0.00, 0.00, 0.00, 0.00},//{0.00, 0.00, 0.75, 0.25}
1109 {0.00, 0.00, 0.00, 0.00},//{0.00, 0.00, 0.50, 0.50}
1110 {0.00, 0.00, 0.00, 0.00},//{0.00, 0.00, 0.25, 0.75}
1111 {0.00, 0.00, 0.00, 1.00}}},
1112 {{{1.00, 0.00, 0.00, 0.00},
1113 {0.75, 0.25, 0.00, 0.00},
1114 {0.50, 0.50, 0.00, 0.00},
1115 {0.25, 0.75, 0.00, 0.00},
1116 {0.00, 1.00, 0.00, 0.00},
1117 {0.75, 0.00, 0.25, 0.00},
1118 {0.75, 0.00, 0.00, 0.25},
1119 {0.50, 0.25, 0.00, 0.25},
1120 {0.25, 0.50, 0.00, 0.25},
1121 {0.00, 0.75, 0.00, 0.25},
1122 {0.50, 0.00, 0.50, 0.00},
1123 {0.50, 0.00, 0.25, 0.25},
1124 {0.50, 0.00, 0.00, 0.50},
1125 {0.25, 0.25, 0.00, 0.50},
1126 {0.00, 0.50, 0.00, 0.50},
1127 {0.25, 0.00, 0.75, 0.00},
1128 {0.25, 0.00, 0.50, 0.25},
1129 {0.25, 0.00, 0.25, 0.50},
1130 {0.25, 0.00, 0.00, 0.75},
1131 {0.00, 0.25, 0.00, 0.75},
1132 {0.00, 0.00, 1.00, 0.00},
1133 {0.00, 0.00, 0.75, 0.25},
1134 {0.00, 0.00, 0.50, 0.50},
1135 {0.00, 0.00, 0.25, 0.75},
1136 {0.00, 0.00, 0.00, 1.00}},
1137 {{1.00, 0.00, 0.00, 0.00},
1138 {0.75, 0.25, 0.00, 0.00},
1139 {0.50, 0.50, 0.00, 0.00},
1140 {0.25, 0.75, 0.00, 0.00},
1141 {0.00, 1.00, 0.00, 0.00},
1142 {1.00, 0.00, 0.00, 0.00},//{0.75, 0.00, 0.25, 0.00}
1143 {0.75, 0.00, 0.00, 0.25},
1144 {0.50, 0.25, 0.00, 0.25},
1145 {0.25, 0.50, 0.00, 0.25},
1146 {0.00, 0.75, 0.00, 0.25},
1147 {1.00, 0.00, 0.00, 0.00},//{0.50, 0.00, 0.50, 0.00}
1148 {0.50/0.75, 0.00, 0.00, 0.25/0.75},//{0.50, 0.00, 0.25, 0.25}
1149 {0.50, 0.00, 0.00, 0.50},
1150 {0.25, 0.25, 0.00, 0.50},
1151 {0.00, 0.50, 0.00, 0.50},
1152 {1.00, 0.00, 0.00, 0.00},//{0.25, 0.00, 0.75, 0.00}
1153 {0.50, 0.00, 0.00, 0.50},//{0.25, 0.00, 0.50, 0.25}
1154 {0.25/0.75, 0.00, 0.00, 0.50/0.75},//{0.25, 0.00, 0.25, 0.50}
1155 {0.25, 0.00, 0.00, 0.75},
1156 {0.00, 0.25, 0.00, 0.75},
1157 {0.00, 0.00, 0.00, 0.00},//{0.00, 0.00, 1.00, 0.00}
1158 {0.00, 0.00, 0.00, 1.00},//{0.00, 0.00, 0.75, 0.25}
1159 {0.00, 0.00, 0.00, 1.00},//{0.00, 0.00, 0.50, 0.50}
1160 {0.00, 0.00, 0.00, 1.00},//{0.00, 0.00, 0.25, 0.75}
1161 {0.00, 0.00, 0.00, 1.00}}}},
1162 {{{{1.00, 0.00, 0.00, 0.00},
1163 {0.75, 0.25, 0.00, 0.00},
1164 {0.50, 0.50, 0.00, 0.00},
1165 {0.25, 0.75, 0.00, 0.00},
1166 {0.00, 1.00, 0.00, 0.00},
1167 {0.75, 0.00, 0.25, 0.00},
1168 {0.75, 0.00, 0.00, 0.25},
1169 {0.50, 0.25, 0.00, 0.25},
1170 {0.25, 0.50, 0.00, 0.25},
1171 {0.00, 0.75, 0.00, 0.25},
1172 {0.50, 0.00, 0.50, 0.00},
1173 {0.50, 0.00, 0.25, 0.25},
1174 {0.50, 0.00, 0.00, 0.50},
1175 {0.25, 0.25, 0.00, 0.50},
1176 {0.00, 0.50, 0.00, 0.50},
1177 {0.25, 0.00, 0.75, 0.00},
1178 {0.25, 0.00, 0.50, 0.25},
1179 {0.25, 0.00, 0.25, 0.50},
1180 {0.25, 0.00, 0.00, 0.75},
1181 {0.00, 0.25, 0.00, 0.75},
1182 {0.00, 0.00, 0.50, 0.50},
1183 {0.00, 0.00, 0.50, 0.50},
1184 {0.00, 0.00, 0.50, 0.50},
1185 {0.00, 0.00, 0.50, 0.50},
1186 {0.00, 0.00, 0.50, 0.50}},
1187 {{1.00, 0.00, 0.00, 0.00},
1188 {0.75, 0.25, 0.00, 0.00},
1189 {0.50, 0.50, 0.00, 0.00},
1190 {0.25, 0.75, 0.00, 0.00},
1191 {0.00, 1.00, 0.00, 0.00},
1192 {0.00, 0.00, 0.00, 0.00},//{0.75, 0.00, 0.25, 0.00}
1193 {0.75, 0.00, 0.00, 0.25},
1194 {0.50, 0.25, 0.00, 0.25},
1195 {0.25, 0.50, 0.00, 0.25},
1196 {0.00, 0.75, 0.00, 0.25},
1197 {0.00, 0.00, 0.00, 0.00},//{0.50, 0.00, 0.50, 0.00}
1198 {0.00, 0.00, 0.00, 0.00},//{0.50, 0.00, 0.25, 0.25}
1199 {0.50, 0.00, 0.00, 0.50},
1200 {0.25, 0.25, 0.00, 0.50},
1201 {0.00, 0.50, 0.00, 0.50},
1202 {0.00, 0.00, 0.00, 0.00},//{0.25, 0.00, 0.75, 0.00}
1203 {0.00, 0.00, 0.00, 0.00},//{0.25, 0.00, 0.50, 0.25}
1204 {0.00, 0.00, 0.00, 0.00},//{0.25, 0.00, 0.25, 0.50}
1205 {0.25, 0.00, 0.00, 0.75},
1206 {0.00, 0.25, 0.00, 0.75},
1207 {0.00, 0.00, 0.00, 0.00},//{0.00, 0.00, 1.00, 0.00}
1208 {0.00, 0.00, 0.00, 0.00},//{0.00, 0.00, 0.75, 0.25}
1209 {0.00, 0.00, 0.00, 0.00},//{0.00, 0.00, 0.50, 0.50}
1210 {0.00, 0.00, 0.00, 0.00},//{0.00, 0.00, 0.25, 0.75}
1211 {0.00, 0.00, 0.00, 0.00}}},//{0.00, 0.00, 0.00, 1.00}
1212 {{{1.00, 0.00, 0.00, 0.00},
1213 {0.75, 0.25, 0.00, 0.00},
1214 {0.50, 0.50, 0.00, 0.00},
1215 {0.25, 0.75, 0.00, 0.00},
1216 {0.00, 1.00, 0.00, 0.00},
1217 {0.75, 0.00, 0.25, 0.00},
1218 {0.75, 0.00, 0.00, 0.25},
1219 {0.50, 0.25, 0.00, 0.25},
1220 {0.25, 0.50, 0.00, 0.25},
1221 {0.00, 0.75, 0.00, 0.25},
1222 {0.50, 0.00, 0.50, 0.00},
1223 {0.50, 0.00, 0.25, 0.25},
1224 {0.50, 0.00, 0.00, 0.50},
1225 {0.25, 0.25, 0.00, 0.50},
1226 {0.00, 0.50, 0.00, 0.50},
1227 {0.25, 0.00, 0.75, 0.00},
1228 {0.25, 0.00, 0.50, 0.25},
1229 {0.25, 0.00, 0.25, 0.50},
1230 {0.25, 0.00, 0.00, 0.75},
1231 {0.00, 0.25, 0.00, 0.75},
1232 {0.00, 0.00, 0.50, 0.50},
1233 {0.00, 0.00, 0.50, 0.50},
1234 {0.00, 0.00, 0.50, 0.50},
1235 {0.00, 0.00, 0.50, 0.50},
1236 {0.00, 0.00, 0.50, 0.50}},
1237 {{1.00, 0.00, 0.00, 0.00},
1238 {0.75, 0.25, 0.00, 0.00},
1239 {0.50, 0.50, 0.00, 0.00},
1240 {0.25, 0.75, 0.00, 0.00},
1241 {0.00, 1.00, 0.00, 0.00},
1242 {1.00, 0.00, 0.00, 0.00},//{0.75, 0.00, 0.25, 0.00}
1243 {0.75, 0.00, 0.00, 0.25},
1244 {0.50, 0.25, 0.00, 0.25},
1245 {0.25, 0.50, 0.00, 0.25},
1246 {0.00, 0.75, 0.00, 0.25},
1247 {1.00, 0.00, 0.00, 0.00},//{0.50, 0.00, 0.50, 0.00}
1248 {0.50/0.75, 0.00, 0.00, 0.25/0.75},//{0.50, 0.00, 0.25, 0.25}
1249 {0.50, 0.00, 0.00, 0.50},
1250 {0.25, 0.25, 0.00, 0.50},
1251 {0.00, 0.50, 0.00, 0.50},
1252 {1.00, 0.00, 0.00, 0.00},//{0.25, 0.00, 0.75, 0.00}
1253 {0.50, 0.00, 0.00, 0.50},//{0.25, 0.00, 0.50, 0.25}
1254 {0.25/0.75, 0.00, 0.00, 0.50/0.75},//{0.25, 0.00, 0.25, 0.50}
1255 {0.25, 0.00, 0.00, 0.75},
1256 {0.00, 0.25, 0.00, 0.75},
1257 {0.00, 0.00, 0.00, 1.00},//{0.00, 0.00, 1.00, 0.00}
1258 {0.00, 0.00, 0.00, 1.00},//{0.00, 0.00, 0.75, 0.25}
1259 {0.00, 0.00, 0.00, 1.00},//{0.00, 0.00, 0.50, 0.50}
1260 {0.00, 0.00, 0.00, 1.00},//{0.00, 0.00, 0.25, 0.75}
1261 {0.00, 0.00, 0.00, 1.00}}}}};
1262
1263 for (size_t j = 0; j < 25; ++j) {
1264 ref_global_tgt_field[j] =
1265 (double)(src_reorder[order_x][order_y][0] + 1) *
1266 ref_weights[at_pole][partial_coverage][with_mask][j][0] +
1267 (double)(src_reorder[order_x][order_y][1] + 1) *
1268 ref_weights[at_pole][partial_coverage][with_mask][j][1] +
1269 (double)(src_reorder[order_x][order_y][2] + 1) *
1270 ref_weights[at_pole][partial_coverage][with_mask][j][2] +
1271 (double)(src_reorder[order_x][order_y][3] + 1) *
1272 ref_weights[at_pole][partial_coverage][with_mask][j][3];
1273 if (ref_global_tgt_field[j] == 0.0)
1274 ref_global_tgt_field[j] = -1.0;
1275 }
1276
1277 if (is_tgt) {
1278 tgt_field = xmalloc(grid_data.num_vertices * sizeof(*tgt_field));
1279 ref_tgt_field = xmalloc(grid_data.num_vertices * sizeof(*ref_tgt_field));
1280 for (size_t i = 0; i < grid_data.num_vertices; ++i)
1281 ref_tgt_field[i] = ref_global_tgt_field[i];
1282 } else {
1283 src_field = xmalloc(grid_data.num_vertices * sizeof(*src_field));
1284 for (size_t i = 0; i < grid_data.num_vertices; ++i)
1285 src_field[i] = (double)(i) + 1.0;
1286 }
1287
1288 yac_interpolation_execute(interpolation, &src_fields, &tgt_field);
1289
1290 // reference data for unstructured grid at the pole is not
1291 // available
1292 if (is_tgt && !(!is_reg2d_grid && at_pole))
1293 for (size_t i = 0; i < grid_data.num_vertices; ++i)
1294 if (fabs(ref_tgt_field[i] - tgt_field[i]) > 1e-3)
1295 PUT_ERR("wrong interpolation result");
1296
1297 free(src_field);
1298 free(tgt_field);
1299 free(ref_tgt_field);
1300 }
1301
1302 yac_interpolation_delete(interpolation);
1303 }
1304
1305 yac_interp_weights_delete(weights);
1306 yac_interp_method_delete(method_stack);
1307 yac_interp_grid_delete(interp_grid);
1308 yac_dist_grid_pair_delete(grid_pair);
1309 yac_basic_grid_delete(grids[1]);
1310 yac_basic_grid_delete(grids[0]);
1311 }
1312 }
1313 }
1314 }
1315 }
1316 }
1317 }
1318
1319 {// Test 8
1320 // barycentric coordinates weighted (with partial coverage) point
1321 // interpolation with two source processes and a single target process
1322
1323 // the global source grid:
1324 // 02------03
1325 // | / |
1326 // | / |
1327 // | / |
1328 // 00------01
1329 //
1330 //---------------
1331 // setup
1332 //---------------
1333
1334 int is_tgt = split_comm_size == 1;
1335 size_t src_num_vertices = 4;
1336 size_t src_num_cells = 2;
1337 int src_num_vertices_per_cell[] = {3, 3};
1338 int src_cell_to_vertex[] = {0,1,3, 0,2,3};
1339 double src_coordinates_x[4] = {-0.5,0.5,-0.5,0.5};
1340 double src_coordinates_y[4] = {-0.5,-0.5,0.5,0.5};
1341 size_t tgt_num_vertices[2] = {5,5};
1342 int tgt_cyclic[2] = {0,0};
1343 double tgt_coordinates_x[5] = {-0.5,-0.25,0.0,0.25,0.5};
1344 double tgt_coordinates_y[5] = {-0.5,-0.25,0.0,0.25,0.5};
1345 int src_cell_mask[4] = {1,1, 0,1};
1346
1347 for (int with_mask = 0; with_mask < 2; ++with_mask) {
1348
1349 struct yac_basic_grid_data grid_data;
1350 if (is_tgt) {
1351 grid_data =
1352 yac_generate_basic_grid_data_reg_2d_deg(
1353 tgt_num_vertices, tgt_cyclic,
1354 tgt_coordinates_x, tgt_coordinates_y);
1355 } else {
1356 grid_data =
1357 yac_generate_basic_grid_data_unstruct_deg(
1358 src_num_vertices, src_num_cells, src_num_vertices_per_cell,
1359 src_coordinates_x, src_coordinates_y, src_cell_to_vertex);
1360 }
1361
1362 struct yac_basic_grid * grids[2] =
1363 {yac_basic_grid_new(grid_names[is_tgt], grid_data),
1364 yac_basic_grid_empty_new(grid_names[is_tgt^1])};
1365
1366 if (with_mask && !is_tgt) {
1367 int * mask = xmalloc(grid_data.num_vertices * sizeof(*mask));
1368 for (size_t i = 0; i < grid_data.num_vertices; ++i)
1369 mask[i] = src_cell_mask[i];
1370 yac_basic_grid_add_mask_nocpy(grids[0], YAC_LOC_CORNER, mask, NULL);
1371 }
1372
1373 struct yac_dist_grid_pair * grid_pair =
1374 yac_dist_grid_pair_new(grids[0], grids[1], MPI_COMM_WORLD);
1375
1376 struct yac_interp_field src_fields[] =
1377 {{.location = YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX,
1378 .masks_idx = (with_mask)?0:SIZE_MAX}};
1379 size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
1380 struct yac_interp_field tgt_field =
1381 {.location = YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX,
1382 .masks_idx = SIZE_MAX};
1383
1384 struct yac_interp_grid * interp_grid =
1385 yac_interp_grid_new(grid_pair, grid_names[0], grid_names[1],
1386 num_src_fields, src_fields, tgt_field);
1387
1388 for (int partial_coverage = 0; partial_coverage < 2; ++partial_coverage) {
1389
1390 struct interp_method * method_stack[] =
1391 {yac_interp_method_avg_new(YAC_INTERP_AVG_BARY, partial_coverage),
1392 yac_interp_method_fixed_new(-1.0), NULL};
1393
1394 struct yac_interp_weights * weights =
1395 yac_interp_method_do_search(method_stack, interp_grid);
1396
1397 for (size_t i = 0; i < num_reorder_types; ++i) {
1398
1399 struct yac_interpolation * interpolation =
1400 yac_interp_weights_get_interpolation(
1401 weights, reorder_types[i], 1,
1402 YAC_FRAC_MASK_NO_VALUE, 1.0, 0.0, NULL, 1, 1);
1403
1404 // check generated interpolation
1405 {
1406 double * src_field = NULL;
1407 double ** src_fields = &src_field;
1408 double * tgt_field = NULL;
1409 double * ref_tgt_field = NULL;
1410 double ref_global_tgt_field[25];
1411 double ref_weights[2][2][25][4] =
1412 {{{{1.00, 0.00, 0.00, 0.00},
1413 {0.75, 0.25, 0.00, 0.00},
1414 {0.50, 0.50, 0.00, 0.00},
1415 {0.25, 0.75, 0.00, 0.00},
1416 {0.00, 1.00, 0.00, 0.00},
1417 {0.75, 0.00, 0.25, 0.00},
1418 {0.75, 0.00, 0.00, 0.25},
1419 {0.50, 0.25, 0.00, 0.25},
1420 {0.25, 0.50, 0.00, 0.25},
1421 {0.00, 0.75, 0.00, 0.25},
1422 {0.50, 0.00, 0.50, 0.00},
1423 {0.50, 0.00, 0.25, 0.25},
1424 {0.50, 0.00, 0.00, 0.50},
1425 {0.25, 0.25, 0.00, 0.50},
1426 {0.00, 0.50, 0.00, 0.50},
1427 {0.25, 0.00, 0.75, 0.00},
1428 {0.25, 0.00, 0.50, 0.25},
1429 {0.25, 0.00, 0.25, 0.50},
1430 {0.25, 0.00, 0.00, 0.75},
1431 {0.00, 0.25, 0.00, 0.75},
1432 {0.00, 0.00, 1.00, 0.00},
1433 {0.00, 0.00, 0.75, 0.25},
1434 {0.00, 0.00, 0.50, 0.50},
1435 {0.00, 0.00, 0.25, 0.75},
1436 {0.00, 0.00, 0.00, 1.00}},
1437 {{1.00, 0.00, 0.00, 0.00},
1438 {0.75, 0.25, 0.00, 0.00},
1439 {0.50, 0.50, 0.00, 0.00},
1440 {0.25, 0.75, 0.00, 0.00},
1441 {0.00, 1.00, 0.00, 0.00},
1442 {0.00, 0.00, 0.00, 0.00},//{0.75, 0.00, 0.25, 0.00}
1443 {0.75, 0.00, 0.00, 0.25},
1444 {0.50, 0.25, 0.00, 0.25},
1445 {0.25, 0.50, 0.00, 0.25},
1446 {0.00, 0.75, 0.00, 0.25},
1447 {0.00, 0.00, 0.00, 0.00},//{0.50, 0.00, 0.50, 0.00}
1448 {0.00, 0.00, 0.00, 0.00},//{0.50, 0.00, 0.25, 0.25}
1449 {0.50, 0.00, 0.00, 0.50},
1450 {0.25, 0.25, 0.00, 0.50},
1451 {0.00, 0.50, 0.00, 0.50},
1452 {0.00, 0.00, 0.00, 0.00},//{0.25, 0.00, 0.75, 0.00}
1453 {0.00, 0.00, 0.00, 0.00},//{0.25, 0.00, 0.50, 0.25}
1454 {0.00, 0.00, 0.00, 0.00},//{0.25, 0.00, 0.25, 0.50}
1455 {0.25, 0.00, 0.00, 0.75},
1456 {0.00, 0.25, 0.00, 0.75},
1457 {0.00, 0.00, 0.00, 0.00},//{0.00, 0.00, 1.00, 0.00}
1458 {0.00, 0.00, 0.00, 0.00},//{0.00, 0.00, 0.75, 0.25}
1459 {0.00, 0.00, 0.00, 0.00},//{0.00, 0.00, 0.50, 0.50}
1460 {0.00, 0.00, 0.00, 0.00},//{0.00, 0.00, 0.25, 0.75}
1461 {0.00, 0.00, 0.00, 1.00}}},
1462 {{{1.00, 0.00, 0.00, 0.00},
1463 {0.75, 0.25, 0.00, 0.00},
1464 {0.50, 0.50, 0.00, 0.00},
1465 {0.25, 0.75, 0.00, 0.00},
1466 {0.00, 1.00, 0.00, 0.00},
1467 {0.75, 0.00, 0.25, 0.00},
1468 {0.75, 0.00, 0.00, 0.25},
1469 {0.50, 0.25, 0.00, 0.25},
1470 {0.25, 0.50, 0.00, 0.25},
1471 {0.00, 0.75, 0.00, 0.25},
1472 {0.50, 0.00, 0.50, 0.00},
1473 {0.50, 0.00, 0.25, 0.25},
1474 {0.50, 0.00, 0.00, 0.50},
1475 {0.25, 0.25, 0.00, 0.50},
1476 {0.00, 0.50, 0.00, 0.50},
1477 {0.25, 0.00, 0.75, 0.00},
1478 {0.25, 0.00, 0.50, 0.25},
1479 {0.25, 0.00, 0.25, 0.50},
1480 {0.25, 0.00, 0.00, 0.75},
1481 {0.00, 0.25, 0.00, 0.75},
1482 {0.00, 0.00, 1.00, 0.00},
1483 {0.00, 0.00, 0.75, 0.25},
1484 {0.00, 0.00, 0.50, 0.50},
1485 {0.00, 0.00, 0.25, 0.75},
1486 {0.00, 0.00, 0.00, 1.00}},
1487 {{1.00, 0.00, 0.00, 0.00},
1488 {0.75, 0.25, 0.00, 0.00},
1489 {0.50, 0.50, 0.00, 0.00},
1490 {0.25, 0.75, 0.00, 0.00},
1491 {0.00, 1.00, 0.00, 0.00},
1492 {1.00, 0.00, 0.00, 0.00},//{0.75, 0.00, 0.25, 0.00}
1493 {0.75, 0.00, 0.00, 0.25},
1494 {0.50, 0.25, 0.00, 0.25},
1495 {0.25, 0.50, 0.00, 0.25},
1496 {0.00, 0.75, 0.00, 0.25},
1497 {1.00, 0.00, 0.00, 0.00},//{0.50, 0.00, 0.50, 0.00}
1498 {0.50/0.75, 0.00, 0.00, 0.25/0.75},//{0.50, 0.00, 0.25, 0.25}
1499 {0.50, 0.00, 0.00, 0.50},
1500 {0.25, 0.25, 0.00, 0.50},
1501 {0.00, 0.50, 0.00, 0.50},
1502 {1.00, 0.00, 0.00, 0.00},//{0.25, 0.00, 0.75, 0.00}
1503 {0.50, 0.00, 0.00, 0.50},//{0.25, 0.00, 0.50, 0.25}
1504 {0.25/0.75, 0.00, 0.00, 0.50/0.75},//{0.25, 0.00, 0.25, 0.50}
1505 {0.25, 0.00, 0.00, 0.75},
1506 {0.00, 0.25, 0.00, 0.75},
1507 {0.00, 0.00, 0.00, 0.00},//{0.00, 0.00, 1.00, 0.00}
1508 {0.00, 0.00, 0.00, 1.00},//{0.00, 0.00, 0.75, 0.25}
1509 {0.00, 0.00, 0.00, 1.00},//{0.00, 0.00, 0.50, 0.50}
1510 {0.00, 0.00, 0.00, 1.00},//{0.00, 0.00, 0.25, 0.75}
1511 {0.00, 0.00, 0.00, 1.00}}}};
1512
1513 for (size_t j = 0; j < 25; ++j) {
1514 ref_global_tgt_field[j] =
1515 1.0 * ref_weights[partial_coverage][with_mask][j][0] +
1516 2.0 * ref_weights[partial_coverage][with_mask][j][1] +
1517 3.0 * ref_weights[partial_coverage][with_mask][j][2] +
1518 4.0 * ref_weights[partial_coverage][with_mask][j][3];
1519 if (ref_global_tgt_field[j] == 0.0)
1520 ref_global_tgt_field[j] = -1.0;
1521 }
1522
1523 if (is_tgt) {
1524 tgt_field = xmalloc(grid_data.num_vertices * sizeof(*tgt_field));
1525 ref_tgt_field = xmalloc(grid_data.num_vertices * sizeof(*ref_tgt_field));
1526 for (size_t i = 0; i < grid_data.num_vertices; ++i)
1527 ref_tgt_field[i] = ref_global_tgt_field[i];
1528 } else {
1529 src_field = xmalloc(grid_data.num_vertices * sizeof(*src_field));
1530 for (size_t i = 0; i < grid_data.num_vertices; ++i)
1531 src_field[i] = (double)(i) + 1.0;
1532 }
1533
1534 yac_interpolation_execute(interpolation, &src_fields, &tgt_field);
1535
1536 if (is_tgt)
1537 for (size_t i = 0; i < grid_data.num_vertices; ++i)
1538 if (fabs(ref_tgt_field[i] - tgt_field[i]) > 1e-3)
1539 PUT_ERR("wrong interpolation result");
1540
1541 free(src_field);
1542 free(tgt_field);
1543 free(ref_tgt_field);
1544 }
1545
1546 yac_interpolation_delete(interpolation);
1547 }
1548
1549 yac_interp_weights_delete(weights);
1550 yac_interp_method_delete(method_stack);
1551 } // partial_coverage
1552 yac_interp_grid_delete(interp_grid);
1553 yac_dist_grid_pair_delete(grid_pair);
1554 yac_basic_grid_delete(grids[1]);
1555 yac_basic_grid_delete(grids[0]);
1556 }
1557 }
1558
1559
1560
1561 {// Test 9
1562 // simple interpolation test with various configuration options
1563
1564 // the global source grid:
1565 // 06--10--07--11--08
1566 // | | |
1567 // 06 02 08 03 09
1568 // | | |
1569 // 03--05--04--07--05
1570 // | | |
1571 // 01 00 03 01 04
1572 // | | |
1573 // 00--00--01--02--02
1574 //
1575 //---------------
1576 // setup
1577 //---------------
1578
1579 int is_tgt = split_comm_size == 1;
1580 double coordinates_x[2][7] = {{0.0,0.1,0.2}, {0.0,0.05,0.1,0.15,0.2}};
1581 double coordinates_y[2][7] = {{0.0,0.1,0.2}, {0.0,0.05,0.1,0.15,0.2}};
1582 size_t const num_cells[2][2] = {{2,2}, {4,4}};
1583 size_t local_start[2][2][2] = {{{0,0},{0,1}}, {{0,0}}};
1584 size_t local_count[2][2][2] = {{{2,1},{2,1}}, {{4,4}}};
1585 int global_mask[2][5*5] = {{1,1,1,
1586 0,1,1,
1587 0,0,1},
1588 {1,1,1,1,1,
1589 1,1,1,1,0,
1590 1,1,1,0,0,
1591 1,1,0,0,0,
1592 1,0,0,0,0}};
1593 int with_halo = 0;
1594 for (size_t i = 0; i <= num_cells[is_tgt][0]; ++i)
1595 coordinates_x[is_tgt][i] *= YAC_RAD;
1596 for (size_t i = 0; i <= num_cells[is_tgt][1]; ++i)
1597 coordinates_y[is_tgt][i] *= YAC_RAD;
1598
1599 struct yac_basic_grid_data grid_data =
1600 utest_generate_basic_grid_data_reg2d(
1601 coordinates_x[is_tgt], coordinates_y[is_tgt], num_cells[is_tgt],
1602 local_start[is_tgt][split_comm_rank],
1603 local_count[is_tgt][split_comm_rank], with_halo);
1604
1605 struct yac_basic_grid * grids[2] =
1606 {yac_basic_grid_new(grid_names[is_tgt], grid_data),
1607 yac_basic_grid_empty_new(grid_names[is_tgt^1])};
1608
1609 int field_mask[5*5];
1610 for (size_t i = 0; i < grid_data.num_vertices; ++i)
1611 field_mask[i] = global_mask[is_tgt][grid_data.vertex_ids[i]];
1612 yac_basic_grid_add_mask(
1613 grids[0], YAC_LOC_CORNER, field_mask, grid_data.num_vertices, NULL);
1614
1615 struct yac_dist_grid_pair * grid_pair =
1616 yac_dist_grid_pair_new(grids[0], grids[1], MPI_COMM_WORLD);
1617
1618 for (int with_src_mask = 0; with_src_mask < 2; ++with_src_mask) {
1619 for (int with_tgt_mask = 0; with_tgt_mask < 2; ++with_tgt_mask) {
1620
1621 struct yac_interp_field src_fields[] =
1622 {{.location = YAC_LOC_CORNER,
1623 .coordinates_idx = SIZE_MAX,
1624 .masks_idx = (with_src_mask)?0:SIZE_MAX}};
1625 size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
1626 struct yac_interp_field tgt_field =
1627 {.location = YAC_LOC_CORNER,
1628 .coordinates_idx = SIZE_MAX,
1629 .masks_idx = (with_tgt_mask)?0:SIZE_MAX};
1630
1631 struct yac_interp_grid * interp_grid =
1632 yac_interp_grid_new(grid_pair, grid_names[0], grid_names[1],
1633 num_src_fields, src_fields, tgt_field);
1634
1635 enum yac_interp_avg_weight_type weight_types[] =
1636 {YAC_INTERP_AVG_ARITHMETIC,
1637 YAC_INTERP_AVG_DIST,
1638 YAC_INTERP_AVG_BARY};
1639 size_t const num_weight_types =
1640 sizeof(weight_types) / sizeof(weight_types[0]);
1641
1642 for (size_t weight_types_idx = 0; weight_types_idx < num_weight_types;
1643 ++weight_types_idx) {
1644 for (int partial_coverage = 0; partial_coverage < 2;
1645 ++partial_coverage) {
1646
1647 struct interp_method * method_stack[] =
1648 {yac_interp_method_avg_new(
1649 weight_types[weight_types_idx], partial_coverage),
1650 yac_interp_method_fixed_new(-1.0), NULL};
1651
1652 struct yac_interp_weights * weights =
1653 yac_interp_method_do_search(method_stack, interp_grid);
1654
1655 for (size_t i = 0; i < num_reorder_types; ++i) {
1656
1657 struct yac_interpolation * interpolation =
1658 yac_interp_weights_get_interpolation(
1659 weights, reorder_types[i], 1,
1660 YAC_FRAC_MASK_NO_VALUE, 1.0, 0.0, NULL, 1, 1);
1661
1662 // check generated interpolation
1663 {
1664 double src_field_data[3*3];
1665 double * src_field = src_field_data;
1666 double ** src_fields = &src_field;
1667 double tgt_field_data[5*5] =
1668 {-2.0,-2.0,-2.0,-2.0,-2.0,
1669 -2.0,-2.0,-2.0,-2.0,-2.0,
1670 -2.0,-2.0,-2.0,-2.0,-2.0,
1671 -2.0,-2.0,-2.0,-2.0,-2.0,
1672 -2.0,-2.0,-2.0,-2.0,-2.0};
1673 double * tgt_field = tgt_field_data;
1674
1675 if (!is_tgt)
1676 for (size_t i = 0; i < grid_data.num_vertices; ++i)
1677 src_field[i] = (double)(grid_data.vertex_ids[i]);
1678
1679 yac_interpolation_execute(
1680 interpolation, &src_fields, &tgt_field);
1681
1682 if (is_tgt) {
1683 size_t src_per_tgt[5*5][4] =
1684 {{0,1,3,4},{0,1,3,4},{0,1,3,4},{1,2,4,5},{1,2,4,5},
1685 {0,1,3,4},{0,1,3,4},{0,1,3,4},{1,2,4,5},{1,2,4,5},
1686 {0,1,3,4},{0,1,3,4},{0,1,3,4},{1,2,4,5},{1,2,4,5},
1687 {3,4,6,7},{3,4,6,7},{3,4,6,7},{4,5,7,8},{4,5,7,8},
1688 {3,4,6,7},{3,4,6,7},{3,4,6,7},{4,5,7,8},{4,5,7,8}};
1689 double A = 0.345491502812526;
1690 double B = 0.154508497187473;
1691 double C = 0.25;
1692 double D = 0.5;
1693 double weights[3][5*5][4] =
1694 {{{C,C,C,C},{C,C,C,C},{C,C,C,C},{C,C,C,C},{C,C,C,C},
1695 {C,C,C,C},{C,C,C,C},{C,C,C,C},{C,C,C,C},{C,C,C,C},
1696 {C,C,C,C},{C,C,C,C},{C,C,C,C},{C,C,C,C},{C,C,C,C},
1697 {C,C,C,C},{C,C,C,C},{C,C,C,C},{C,C,C,C},{C,C,C,C},
1698 {C,C,C,C},{C,C,C,C},{C,C,C,C},{C,C,C,C},{C,C,C,C}},
1699 {{1,0,0,0},{A,A,B,B},{0,1,0,0},{A,A,B,B},{0,1,0,0},
1700 {A,B,A,B},{C,C,C,C},{B,A,B,A},{C,C,C,C},{B,A,B,A},
1701 {0,0,1,0},{B,B,A,A},{0,0,0,1},{B,B,A,A},{0,0,0,1},
1702 {A,B,A,B},{C,C,C,C},{B,A,B,A},{C,C,C,C},{B,A,B,A},
1703 {0,0,1,0},{B,B,A,A},{0,0,0,1},{B,B,A,A},{0,0,0,1}},
1704 {{1,0,0,0},{D,D,0,0},{0,1,0,0},{D,D,0,0},{0,1,0,0},
1705 {D,0,D,0},{D,0,0,D},{0,D,0,D},{D,0,0,D},{0,D,0,D},
1706 {0,0,1,0},{0,0,D,D},{0,0,0,1},{0,0,D,D},{0,0,0,1},
1707 {D,0,D,0},{D,0,0,D},{0,D,0,D},{D,0,0,D},{0,D,0,D},
1708 {0,0,1,0},{0,0,D,D},{0,0,0,1},{0,0,D,D},{0,0,0,1}}};
1709
1710 for (size_t i = 0; i < grid_data.num_vertices; ++i) {
1711
1712 size_t tgt_idx = grid_data.vertex_ids[i];
1713 double src_values[4];
1714 double w[4];
1715
1716 double ref_tgt_value = 0.0;
1717
1718 int contains_masked_src = 0;
1719 int contains_unmasked_src = 0;
1720
1721 for (int j = 0; j < 4; ++j) {
1722 size_t src_idx = src_per_tgt[tgt_idx][j];
1723 w[j] = weights[weight_types_idx][tgt_idx][j];
1724 if (with_src_mask &&
1725 !global_mask[0][src_idx] &&
1726 (w[j] != 0.0)) {
1727 src_values[j] = -1.0;
1728 w[j] = 0.0;
1729 contains_masked_src = 1;
1730 } else {
1731 src_values[j] = (double)src_idx;
1732 contains_unmasked_src |= w[j] != 0.0;
1733 }
1734 }
1735
1736 if (with_tgt_mask && !global_mask[1][tgt_idx]) {
1737 ref_tgt_value = -2.0;
1738 } else if (contains_masked_src) {
1739 if (partial_coverage && contains_unmasked_src) {
1740 yac_correct_weights(4, w);
1741 for (int j = 0; j < 4; ++j)
1742 ref_tgt_value += w[j] * src_values[j];
1743 } else {
1744 ref_tgt_value = -1.0;
1745 }
1746 } else {
1747 for (int j = 0; j < 4; ++j)
1748 ref_tgt_value += w[j] * src_values[j];
1749 }
1750
1751 // some special cases
1752 if (with_src_mask && partial_coverage &&
1753 (weight_types[weight_types_idx] ==
1754 YAC_INTERP_AVG_DIST) &&
1755 (!with_tgt_mask || global_mask[1][tgt_idx])) {
1756 switch (tgt_idx) {
1757 case(10):
1758 ref_tgt_value = 1.73879612503;
1759 break;
1760 case(20):
1761 case(22):
1762 ref_tgt_value = 4.0;
1763 default:
1764 break;
1765 }
1766 }
1767
1768 if (fabs(ref_tgt_value - tgt_field[i]) > 1e-3)
1769 PUT_ERR("wrong interpolation result");
1770 }
1771 }
1772 }
1773
1774 yac_interpolation_delete(interpolation);
1775 } // reorder type
1776 yac_interp_weights_delete(weights);
1777 yac_interp_method_delete(method_stack);
1778 } // partial coverage
1779 } // weight type
1780 yac_interp_grid_delete(interp_grid);
1781 } // with tgt mask
1782 } // with src mask
1783 yac_dist_grid_pair_delete(grid_pair);
1784 yac_basic_grid_delete(grids[1]);
1785 yac_basic_grid_delete(grids[0]);
1786 }
1787
1788 {// Test 10
1789 // barycentric coordinates weighted (cell touches pole and has duplicated
1790 // longitude coordinate at pole)
1791 // interpolation with two source processes and a single target process
1792
1793 // the global source grid is a 1x1 grid:
1794 // 02--03--03
1795 // | |
1796 // 01 00 02
1797 // | |
1798 // 00--00--01
1799 //
1800 //---------------
1801 // setup
1802 //---------------
1803
1804 int is_tgt = split_comm_size == 1;
1805 double src_coordinates_x[2][2] = {{-0.5,0.5}, {0.5,-0.5}};
1806 double src_coordinates_y[2][2] = {{89,90}, {90,89}};
1807 size_t const src_num_cells[2] = {1,1};
1808 size_t src_local_start[2] = {0,0};
1809 size_t src_local_count[2] = {1,1};
1810 yac_int src_reorder[2][2][4] =
1811 {{{0,1,2,3},{2,3,0,1}},{{1,0,3,2},{3,2,1,0}}};
1812 double tgt_coordinates_x[5] = {-0.5,-0.25,0.0,0.25,0.5};
1813 double tgt_coordinates_y[5] = {89.00,89.25,89.50,89.75,90.00};
1814 size_t const tgt_num_cells[2] = {4,4};
1815 size_t tgt_local_start[2] = {0,0};
1816 size_t tgt_local_count[2] = {4,4};
1817
1818 for (size_t i = 0; i < 2 * 2; ++i)
1819 (&(src_coordinates_x[0][0]))[i] *= YAC_RAD;
1820 for (size_t i = 0; i < 2 * 2; ++i)
1821 (&(src_coordinates_y[0][0]))[i] *= YAC_RAD;
1822 for (size_t i = 0; i < 5; ++i) tgt_coordinates_x[i] *= YAC_RAD;
1823 for (size_t i = 0; i < 5; ++i)
1824 (&(tgt_coordinates_y[0]))[i] *= YAC_RAD;
1825
1826 for (int order_x = 0; order_x < 2; ++order_x) {
1827 for (int order_y = 0; order_y < 2; ++order_y) {
1828
1829 struct yac_basic_grid_data grid_data;
1830
1831 if (is_tgt) {
1832 grid_data =
1833 utest_generate_basic_grid_data_reg2d(
1834 tgt_coordinates_x, tgt_coordinates_y,
1835 tgt_num_cells, tgt_local_start, tgt_local_count, 0);
1836 } else {
1837 grid_data =
1838 utest_generate_basic_grid_data_reg2d(
1839 src_coordinates_x[order_x],
1840 src_coordinates_y[order_y], src_num_cells,
1841 src_local_start, src_local_count, 0);
1842 }
1843
1844 if (!is_tgt) {
1845
1846 yac_int * vertex_ids =
1847 xmalloc(grid_data.num_vertices * sizeof(*vertex_ids));
1848 for (size_t i = 0; i < grid_data.num_vertices; ++i)
1849 vertex_ids[i] = (yac_int)(src_reorder[order_x][order_y][i]);
1850 grid_data.vertex_ids = vertex_ids;
1851 }
1852
1853 struct yac_basic_grid * grids[2] =
1854 {yac_basic_grid_new(grid_names[is_tgt], grid_data),
1855 yac_basic_grid_empty_new(grid_names[is_tgt^1])};
1856
1857 if (!is_tgt) {
1858
1859 yac_coordinate_pointer vertex_coordinates =
1860 xmalloc(grid_data.num_vertices * sizeof(*vertex_coordinates));
1861 memcpy(
1862 vertex_coordinates, grid_data.vertex_coordinates,
1863 grid_data.num_vertices * sizeof(*vertex_coordinates));
1864 // set x and y coordinate of all pole points to 0.0 in order to
1865 // reproduce an error that occured in ICON
1866 for (size_t i = 0; i < grid_data.num_vertices; ++i) {
1867 if (fabs(1.0 - fabs(vertex_coordinates[i][2])) < 1e-9) {
1868 vertex_coordinates[i][0] = 0.0;
1869 vertex_coordinates[i][1] = 0.0;
1870 }
1871 }
1872 yac_basic_grid_add_coordinates_nocpy(
1873 grids[0], YAC_LOC_CORNER, vertex_coordinates);
1874 }
1875
1876 struct yac_dist_grid_pair * grid_pair =
1877 yac_dist_grid_pair_new(grids[0], grids[1], MPI_COMM_WORLD);
1878
1879 struct yac_interp_field src_fields[] =
1880 {{.location = YAC_LOC_CORNER, .coordinates_idx = 0,
1881 .masks_idx = SIZE_MAX}};
1882 size_t num_src_fields =
1883 sizeof(src_fields) / sizeof(src_fields[0]);
1884 struct yac_interp_field tgt_field =
1885 {.location = YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX,
1886 .masks_idx = SIZE_MAX};
1887
1888 struct yac_interp_grid * interp_grid =
1889 yac_interp_grid_new(grid_pair, grid_names[0], grid_names[1],
1890 num_src_fields, src_fields, tgt_field);
1891
1892 int partial_coverage = 0;
1893 struct interp_method * method_stack[] =
1894 {yac_interp_method_avg_new(
1895 YAC_INTERP_AVG_BARY, partial_coverage),
1896 yac_interp_method_fixed_new(-1.0), NULL};
1897
1898 struct yac_interp_weights * weights =
1899 yac_interp_method_do_search(method_stack, interp_grid);
1900
1901 for (size_t i = 0; i < num_reorder_types; ++i) {
1902
1903 struct yac_interpolation * interpolation =
1904 yac_interp_weights_get_interpolation(
1905 weights, reorder_types[i], 1,
1906 YAC_FRAC_MASK_NO_VALUE, 1.0, 0.0, NULL, 1, 1);
1907
1908 // check generated interpolation
1909 {
1910 double * src_field = NULL;
1911 double ** src_fields = &src_field;
1912 double * tgt_field = NULL;
1913 double * ref_tgt_field = NULL;
1914 double ref_global_tgt_field[25];
1915 double ref_weights[25][4] =
1916 {{1.00, 0.00, 0.00, 0.00},
1917 {0.75, 0.25, 0.00, 0.00},
1918 {0.50, 0.50, 0.00, 0.00},
1919 {0.25, 0.75, 0.00, 0.00},
1920 {0.00, 1.00, 0.00, 0.00},
1921 {0.75, 0.00, 0.25, 0.00},
1922 {0.75, 0.00, 0.00, 0.25},
1923 {0.50, 0.25, 0.00, 0.25},
1924 {0.25, 0.50, 0.00, 0.25},
1925 {0.00, 0.75, 0.00, 0.25},
1926 {0.50, 0.00, 0.50, 0.00},
1927 {0.50, 0.00, 0.25, 0.25},
1928 {0.50, 0.00, 0.00, 0.50},
1929 {0.25, 0.25, 0.00, 0.50},
1930 {0.00, 0.50, 0.00, 0.50},
1931 {0.25, 0.00, 0.75, 0.00},
1932 {0.25, 0.00, 0.50, 0.25},
1933 {0.25, 0.00, 0.25, 0.50},
1934 {0.25, 0.00, 0.00, 0.75},
1935 {0.00, 0.25, 0.00, 0.75},
1936 {0.00, 0.00, 0.50, 0.50},
1937 {0.00, 0.00, 0.50, 0.50},
1938 {0.00, 0.00, 0.50, 0.50},
1939 {0.00, 0.00, 0.50, 0.50},
1940 {0.00, 0.00, 0.50, 0.50}};
1941
1942 for (size_t j = 0; j < 25; ++j) {
1943 ref_global_tgt_field[j] =
1944 (double)(src_reorder[order_x][order_y][0] + 1) *
1945 ref_weights[j][0] +
1946 (double)(src_reorder[order_x][order_y][1] + 1) *
1947 ref_weights[j][1] +
1948 (double)(src_reorder[order_x][order_y][2] + 1) *
1949 ref_weights[j][2] +
1950 (double)(src_reorder[order_x][order_y][3] + 1) *
1951 ref_weights[j][3];
1952 if (ref_global_tgt_field[j] == 0.0)
1953 ref_global_tgt_field[j] = -1.0;
1954 }
1955
1956 if (is_tgt) {
1957 tgt_field = xmalloc(grid_data.num_vertices * sizeof(*tgt_field));
1958 ref_tgt_field = xmalloc(grid_data.num_vertices * sizeof(*ref_tgt_field));
1959 for (size_t i = 0; i < grid_data.num_vertices; ++i)
1960 ref_tgt_field[i] = ref_global_tgt_field[i];
1961 } else {
1962 src_field = xmalloc(grid_data.num_vertices * sizeof(*src_field));
1963 for (size_t i = 0; i < grid_data.num_vertices; ++i)
1964 src_field[i] = (double)(i) + 1.0;
1965 }
1966
1967 yac_interpolation_execute(interpolation, &src_fields, &tgt_field);
1968
1969 if (is_tgt)
1970 for (size_t i = 0; i < grid_data.num_vertices; ++i)
1971 if (fabs(ref_tgt_field[i] - tgt_field[i]) > 1e-3)
1972 PUT_ERR("wrong interpolation result");
1973
1974 free(src_field);
1975 free(tgt_field);
1976 free(ref_tgt_field);
1977 }
1978
1979 yac_interpolation_delete(interpolation);
1980 }
1981
1982 yac_interp_weights_delete(weights);
1983 yac_interp_method_delete(method_stack);
1984 yac_interp_grid_delete(interp_grid);
1985 yac_dist_grid_pair_delete(grid_pair);
1986 yac_basic_grid_delete(grids[1]);
1987 yac_basic_grid_delete(grids[0]);
1988 }
1989 }
1990 }
1991
1992 MPI_Comm_free(&split_comm);
1993
1994 xt_finalize();
1995
1996 MPI_Finalize();
1997
1998 return TEST_EXIT_CODE;
1999}
yac_basic_grid_new
struct yac_basic_grid * yac_basic_grid_new(char const *name, struct yac_basic_grid_data grid_data)
Definition basic_grid.c:57
yac_basic_grid_add_mask
size_t yac_basic_grid_add_mask(struct yac_basic_grid *grid, enum yac_location location, int const *mask, size_t count, char const *mask_name)
Definition basic_grid.c:266
yac_basic_grid_add_coordinates_nocpy
size_t yac_basic_grid_add_coordinates_nocpy(struct yac_basic_grid *grid, enum yac_location location, yac_coordinate_pointer coordinates)
Definition basic_grid.c:202
yac_basic_grid_empty_new
struct yac_basic_grid * yac_basic_grid_empty_new(char const *name)
Definition basic_grid.c:70
yac_basic_grid_delete
void yac_basic_grid_delete(struct yac_basic_grid *grid)
Definition basic_grid.c:77
yac_basic_grid_add_mask_nocpy
size_t yac_basic_grid_add_mask_nocpy(struct yac_basic_grid *grid, enum yac_location location, int const *mask, char const *mask_name)
Definition basic_grid.c:244
yac_generate_basic_grid_data_unstruct_deg
struct yac_basic_grid_data yac_generate_basic_grid_data_unstruct_deg(size_t nbr_vertices, size_t nbr_cells, int *num_vertices_per_cell, double *x_vertices, double *y_vertices, int *cell_to_vertex)
Definition grid_unstruct.c:270
yac_generate_basic_grid_data_unstruct
struct yac_basic_grid_data yac_generate_basic_grid_data_unstruct(size_t nbr_vertices, size_t nbr_cells, int *num_vertices_per_cell, double *x_vertices, double *y_vertices, int *cell_to_vertex)
Definition grid_unstruct.c:260
yac_generate_basic_grid_data_reg_2d_deg
struct yac_basic_grid_data yac_generate_basic_grid_data_reg_2d_deg(size_t nbr_vertices[2], int cyclic[2], double *lon_vertices, double *lat_vertices)
Definition grid_reg2d.c:71
yac_correct_weights
void yac_correct_weights(size_t nSourceCells, double *weight)
correct interpolation weights
Definition clipping.c:1513
yac_dist_grid_pair_delete
void yac_dist_grid_pair_delete(struct yac_dist_grid_pair *grid_pair)
Definition dist_grid.c:2377
yac_dist_grid_pair_new
struct yac_dist_grid_pair * yac_dist_grid_pair_new(struct yac_basic_grid *grid_a, struct yac_basic_grid *grid_b, MPI_Comm comm)
Definition dist_grid.c:2089
dist_grid_utils.h
YAC_RAD
#define YAC_RAD
Definition dynamic_config.c:12
LLtoXYZ_deg
static void LLtoXYZ_deg(double lon, double lat, double p_out[])
Definition geometry.h:278
normalise_vector
static void normalise_vector(double v[])
Definition geometry.h:754
yac_interp_grid_delete
void yac_interp_grid_delete(struct yac_interp_grid *interp_grid)
Definition interp_grid.c:556
yac_interp_grid_new
struct yac_interp_grid * yac_interp_grid_new(struct yac_dist_grid_pair *grid_pair, char const *src_grid_name, char const *tgt_grid_name, size_t num_src_fields, struct yac_interp_field const *src_fields, struct yac_interp_field const tgt_field)
Definition interp_grid.c:31
yac_interp_method_delete
void yac_interp_method_delete(struct interp_method **method)
Delete an interpolation stack and free its resources (but not the pointer array).
Definition interp_method.c:62
yac_interp_method_do_search
struct yac_interp_weights * yac_interp_method_do_search(struct interp_method **method, struct yac_interp_grid *interp_grid)
Perform weight computation using given interpolation stack and grid.
Definition interp_method.c:16
interp_method.h
Defines the interface of the interpolation method "base class" in YAC.
yac_interp_method_avg_new
struct interp_method * yac_interp_method_avg_new(enum yac_interp_avg_weight_type weight_type, int partial_coverage)
Definition interp_method_avg.c:1005
interp_method_avg.h
yac_interp_avg_weight_type
yac_interp_avg_weight_type
Definition interp_method_avg.h:12
YAC_INTERP_AVG_DIST
@ YAC_INTERP_AVG_DIST
Definition interp_method_avg.h:14
YAC_INTERP_AVG_ARITHMETIC
@ YAC_INTERP_AVG_ARITHMETIC
Definition interp_method_avg.h:13
YAC_INTERP_AVG_BARY
@ YAC_INTERP_AVG_BARY
Definition interp_method_avg.h:15
yac_interp_method_fixed_new
struct interp_method * yac_interp_method_fixed_new(double value)
Definition interp_method_fixed.c:63
interp_method_fixed.h
yac_interp_weights_get_interpolation
struct yac_interpolation * yac_interp_weights_get_interpolation(struct yac_interp_weights *weights, enum yac_interp_weights_reorder_type reorder, size_t collection_size, double frac_mask_fallback_value, double scaling_factor, double scaling_summand, char const *yaxt_exchanger_name, int is_source, int is_target)
Definition interp_weights.c:4885
yac_interp_weights_delete
void yac_interp_weights_delete(struct yac_interp_weights *weights)
Definition interp_weights.c:6459
yac_interp_weights_reorder_type
yac_interp_weights_reorder_type
Definition interp_weights.h:44
YAC_MAPPING_ON_TGT
@ YAC_MAPPING_ON_TGT
weights will be applied at target processes
Definition interp_weights.h:46
YAC_MAPPING_ON_SRC
@ YAC_MAPPING_ON_SRC
weights will be applied at source processes
Definition interp_weights.h:45
yac_interpolation_execute
void yac_interpolation_execute(struct yac_interpolation *interp, double ***src_fields, double **tgt_field)
Execute interpolation synchronously and write results to the target field.
Definition interpolation.c:247
yac_interpolation_delete
void yac_interpolation_delete(struct yac_interpolation *interp)
Free an interpolation object and release all resources.
Definition interpolation.c:364
YAC_FRAC_MASK_NO_VALUE
double const YAC_FRAC_MASK_NO_VALUE
Definition interpolation.c:27
YAC_LOC_CORNER
@ YAC_LOC_CORNER
Definition location.h:15
YAC_LOC_CELL
@ YAC_LOC_CELL
Definition location.h:14
test_core.weights
weights
Definition test_core.py:69
test_core.src_field
src_field
Definition test_core.py:61
test_raw_exchange.i
int i
Definition test_raw_exchange.py:89
test_raw_exchange.tgt_field
tgt_field
Definition test_raw_exchange.py:64
yac_generate_weights.src_mask
src_mask
Definition yac_generate_weights.py:167
xmalloc
#define xmalloc(size)
Definition ppm_xfuncs.h:66
grid_data
Definition test_dummy_coupling3_c.c:36
grid_data::num_vertices_per_cell
int * num_vertices_per_cell
Definition test_dummy_coupling3_c.c:42
grid_data::cell_to_vertex
int * cell_to_vertex
Definition test_dummy_coupling3_c.c:43
interp_method
Definition interp_method_internal.h:36
yac_basic_grid_data
Definition basic_grid_data.h:13
yac_basic_grid_data::vertex_ids
yac_int * vertex_ids
Definition basic_grid_data.h:16
yac_basic_grid
Definition basic_grid.c:24
yac_dist_grid_pair
Definition dist_grid.c:128
yac_interp_field
Definition basic_grid.h:15
yac_interp_field::location
enum yac_location location
Definition basic_grid.h:16
yac_interp_grid
Definition interp_grid.c:22
yac_interp_grid::tgt_field
struct yac_interp_field tgt_field
Definition interp_grid.c:26
yac_interp_grid::num_src_fields
size_t num_src_fields
Definition interp_grid.c:27
yac_interp_grid::grid_pair
struct yac_dist_grid_pair * grid_pair
Definition interp_grid.c:25
yac_interp_grid::src_fields
struct yac_interp_field src_fields[]
Definition interp_grid.c:28
yac_interp_weights
Definition interp_weights.c:122
yac_interpolation
Definition interpolation.c:30
test_common.h
grid_names
static char const * grid_names[2]
Definition test_interp_method_avg_parallel.c:29
num_reorder_types
size_t num_reorder_types
Definition test_interp_method_avg_parallel.c:28
reorder_types
enum yac_interp_weights_reorder_type reorder_types[]
Definition test_interp_method_avg_parallel.c:26
split_comm
static MPI_Comm split_comm
Definition test_interp_method_file_parallel.c:55
split_comm_rank
int split_comm_rank
Definition test_interp_method_fixed_parallel.c:28
split_comm_size
int split_comm_size
Definition test_interp_method_fixed_parallel.c:28
num_weight_types
size_t num_weight_types
Definition test_interp_method_spmap_parallel.c:32
weight_types
enum yac_interp_spmap_weight_type weight_types[]
Definition test_interp_method_spmap_parallel.c:28
tgt_field_data
static double tgt_field_data[MAX_COLLECTION_SIZE][NUM_TGT_POINTS]
Definition test_interp_operator_fixed.c:31
coordinates_x
double coordinates_x[]
Definition test_interpolation_parallel3.c:105
num_cells
size_t num_cells[2]
Definition test_interpolation_parallel3.c:104
cell_coordinates_y
double cell_coordinates_y[]
Definition test_interpolation_parallel3.c:108
coordinates_y
double coordinates_y[]
Definition test_interpolation_parallel3.c:106
ref_weights
double ref_weights[4 *36]
Definition test_interpolation_parallel3.c:64
with_halo
int with_halo
Definition test_interpolation_parallel3.c:110
cell_coordinates_x
double cell_coordinates_x[]
Definition test_interpolation_parallel3.c:107
mask
static int mask[16]
Definition test_read_icon_common.h:34
TEST_EXIT_CODE
#define TEST_EXIT_CODE
Definition tests.h:15
PUT_ERR
#define PUT_ERR(string)
Definition tests.h:10
yac_mpi.h
yac_int
YAC_INT yac_int
Definition yac_types.h:15
yac_coordinate_pointer
double(* yac_coordinate_pointer)[3]
Definition yac_types.h:21