YAC 3.21.0
Yet Another Coupler
Loading...
Searching...
No Matches
test_interp_method_ncc_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
7#include "tests.h"
8#include "test_common.h"
9#include "geometry.h"
10#include "interpolation/methods/interp_method.h"
11#include "interpolation/methods/interp_method_ncc.h"
12#include "interpolation/methods/interp_method_fixed.h"
13#include "dist_grid_utils.h"
14#include "clipping.h"
15#include "yac_mpi.h"
16
17#include <mpi.h>
18#include <yaxt.h>
19
27size_t num_reorder_types = sizeof(reorder_types) / sizeof(reorder_types[0]);
28static char const * grid_names[2] = {"src_grid", "tgt_grid"};
29
30static double utest_compute_ref_tgt_value(
31 size_t tgt_corner, int * tgt_global_core_mask, int * tgt_global_field_mask,
32 int * src_global_core_mask, int * src_global_field_maske, int with_src_mask,
33 int with_tgt_mask, enum yac_interp_ncc_weight_type weight_type,
34 int partial_coverage);
35
36int main(void) {
37
38 MPI_Init(NULL, NULL);
39
40 xt_initialize(MPI_COMM_WORLD);
41
42 int comm_rank, comm_size;
43 MPI_Comm_rank(MPI_COMM_WORLD, &comm_rank);
44 MPI_Comm_size(MPI_COMM_WORLD, &comm_size);
45 MPI_Barrier(MPI_COMM_WORLD);
46
47 if (comm_size != 3) {
48 PUT_ERR("ERROR: wrong number of processes");
49 xt_finalize();
50 MPI_Finalize();
51 return TEST_EXIT_CODE;
52 }
53
54 MPI_Comm split_comm;
55 MPI_Comm_split(
56 MPI_COMM_WORLD, comm_rank < 2, 0, &split_comm);
57
58 int split_comm_rank, split_comm_size;
59 MPI_Comm_rank(split_comm, &split_comm_rank);
60 MPI_Comm_size(split_comm, &split_comm_size);
61
62 {// test with various configuration options
63
64 // The global source grid has 5x4 cells:
65 //
66 // 24--44--25--45--26--46--27--47--28--48--29
67 // | | | | | |
68 // 34 15 36 16 38 17 40 18 42 19 43
69 // | | | | | |
70 // 18--33--19--35--20--37--21--39--22--41--23
71 // | | | | | |
72 // 23 10 25 11 27 12 29 13 31 14 32
73 // | | | | | |
74 // 12--22--13--24--14--26--15--28--16--30--17
75 // | | | | | |
76 // 12 05 14 06 16 07 18 08 20 09 21
77 // | | | | | |
78 // 06--11--07--13--08--15--09--17--10--19--11
79 // | | | | | |
80 // 01 00 03 01 05 02 07 03 09 04 10
81 // | | | | | |
82 // 00--01--01--02--02--04--03--06--04--08--05
83 //
84 //---------------
85 // setup
86 //---------------
87
88 int is_tgt = split_comm_size == 1;
89 double coordinates_x[2][10] =
90 {{0,1,2,3,4,5}, {0.25,0.75,1.25,1.75,2.25,2.75,3.25,3.75,4.25,4.75}};
91 double coordinates_y[2][8] =
92 {{-2,-1,0,1,2}, {-1.75,-1.25,-0.75,-0.25,0.25,0.75,1.25,1.75}};
93 double cell_coordiantes_x[5 * 4] =
94 {0.5,1.5,2.5,3.5,4.5,
95 0.5,1.5,2.5,3.5,4.5,
96 0.5,1.5,2.5,3.5,4.5,
97 0.5,1.5,2.5,3.5,4.5};
98 double cell_coordiantes_y[5 * 4] =
99 {-1.5,-1.5,-1.5,-1.5,-1.5,
100 -0.5,-0.5,-0.5,-0.5,-0.5,
101 0.5,0.5,0.5,0.5,0.5,
102 1.5,1.5,1.5,1.5,1.5};
103 size_t const num_cells[2][2] = {{5,4}, {9,7}};
104 size_t local_start[2][2][2] = {{{0,0},{3,0}}, {{0,0}}};
105 size_t local_count[2][2][2] = {{{3,4},{2,4}}, {{9,7}}};
106 int global_core_mask[2][10*8] =
107 {{1,0,0,0,0,
108 1,1,0,0,0,
109 1,1,1,1,0,
110 1,1,1,1,0},
111 {1,1,1,0,0,0,0,0,0,0,
112 1,1,1,0,0,0,0,0,0,0,
113 1,1,1,1,1,1,0,0,0,0,
114 1,1,1,1,1,1,0,0,0,0,
115 1,1,1,1,1,1,1,1,1,0,
116 1,1,1,1,1,1,1,1,1,0,
117 1,1,1,1,1,1,1,1,1,0,
118 1,1,1,1,1,1,1,1,1,0}};
119 int global_field_mask[2][10*8] =
120 {{0,1,1,1,1,
121 0,1,1,1,1,
122 0,1,1,1,1,
123 0,1,1,1,1},
124 {1,1,1,1,1,1,1,1,1,1,
125 1,1,1,1,1,1,1,1,1,1,
126 1,1,1,1,1,1,1,1,1,1,
127 1,1,1,1,1,1,1,1,1,1,
128 1,1,1,1,1,1,1,1,1,1,
129 1,1,1,1,1,1,1,1,1,1,
130 1,1,1,1,1,1,1,1,1,1,
131 0,0,0,0,0,0,0,0,0,0}};
132 int with_halo = 0;
133 for (size_t i = 0; i <= num_cells[is_tgt][0]; ++i)
134 coordinates_x[is_tgt][i] *= YAC_RAD;
135 for (size_t i = 0; i <= num_cells[is_tgt][1]; ++i)
136 coordinates_y[is_tgt][i] *= YAC_RAD;
137
138 struct yac_basic_grid_data grid_data =
139 utest_generate_basic_grid_data_reg2d(
140 coordinates_x[is_tgt], coordinates_y[is_tgt], num_cells[is_tgt],
141 local_start[is_tgt][split_comm_rank],
142 local_count[is_tgt][split_comm_rank], with_halo);
143
144 if (is_tgt) {
145 for (size_t i = 0; i < grid_data.num_vertices; ++i)
146 grid_data.core_vertex_mask[i] =
147 global_core_mask[is_tgt][grid_data.vertex_ids[i]];
148 free(grid_data.core_edge_mask);
149 free(grid_data.core_cell_mask);
150 grid_data.core_edge_mask = NULL;
151 grid_data.core_cell_mask = NULL;
152 } else {
153 for (size_t i = 0; i < grid_data.num_cells; ++i)
154 grid_data.core_cell_mask[i] =
155 global_core_mask[is_tgt][grid_data.cell_ids[i]];
156 }
157
158 struct yac_basic_grid * grids[2] =
159 {yac_basic_grid_new(grid_names[is_tgt], grid_data),
160 yac_basic_grid_empty_new(grid_names[is_tgt^1])};
161
162 int field_mask[10*8];
163 double field_coords[10*8][3];
164 if (is_tgt) {
165 for (size_t i = 0; i < grid_data.num_vertices; ++i)
166 field_mask[i] = global_field_mask[is_tgt][grid_data.vertex_ids[i]];
167 yac_basic_grid_add_mask(
168 grids[0], YAC_LOC_CORNER, field_mask, grid_data.num_vertices, NULL);
169 } else {
170 for (size_t i = 0; i < grid_data.num_cells; ++i) {
171 field_mask[i] = global_field_mask[is_tgt][grid_data.cell_ids[i]];
172 LLtoXYZ_deg(
173 cell_coordiantes_x[grid_data.cell_ids[i]],
174 cell_coordiantes_y[grid_data.cell_ids[i]], field_coords[i]);
175 }
176 yac_basic_grid_add_mask(
177 grids[0], YAC_LOC_CELL, field_mask, grid_data.num_cells, NULL);
178 yac_basic_grid_add_coordinates(
179 grids[0], YAC_LOC_CELL, field_coords, grid_data.num_cells);
180 }
181
182 struct yac_dist_grid_pair * grid_pair =
183 yac_dist_grid_pair_new(grids[0], grids[1], MPI_COMM_WORLD);
184
185 for (int with_src_mask = 0; with_src_mask < 2; ++with_src_mask) {
186 for (int with_tgt_mask = 0; with_tgt_mask < 2; ++with_tgt_mask) {
187
188 struct yac_interp_field src_fields[] =
189 {{.location = YAC_LOC_CELL,
190 .coordinates_idx = 0,
191 .masks_idx = (with_src_mask)?0:SIZE_MAX}};
192 size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
193 struct yac_interp_field tgt_field =
194 {.location = YAC_LOC_CORNER,
195 .coordinates_idx = SIZE_MAX,
196 .masks_idx = (with_tgt_mask)?0:SIZE_MAX};
197
198 struct yac_interp_grid * interp_grid =
199 yac_interp_grid_new(grid_pair, grid_names[0], grid_names[1],
200 num_src_fields, src_fields, tgt_field);
201
202 enum yac_interp_ncc_weight_type weight_types[] =
203 {YAC_INTERP_NCC_AVG,
204 YAC_INTERP_NCC_DIST};
205 size_t const num_weight_types =
206 sizeof(weight_types) / sizeof(weight_types[0]);
207
208 for (size_t weight_types_idx = 0; weight_types_idx < num_weight_types;
209 ++weight_types_idx) {
210 for (int partial_coverage = 0; partial_coverage < 2;
211 ++partial_coverage) {
212
213 struct interp_method * method_stack[] =
214 {yac_interp_method_ncc_new(
215 weight_types[weight_types_idx], partial_coverage),
216 yac_interp_method_fixed_new(-1.0), NULL};
217
218 struct yac_interp_weights * weights =
219 yac_interp_method_do_search(method_stack, interp_grid);
220
221 for (size_t i = 0; i < num_reorder_types; ++i) {
222
223 struct yac_interpolation * interpolation =
224 yac_interp_weights_get_interpolation(
225 weights, reorder_types[i], 1,
226 YAC_FRAC_MASK_NO_VALUE, 1.0, 0.0, NULL, 1, 1);
227
228 // check generated interpolation
229 {
230 double src_field_data[5*4];
231 double * src_field = src_field_data;
232 double ** src_fields = &src_field;
233 double tgt_field_data[10*8] =
234 {-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,
235 -2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,
236 -2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,
237 -2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,
238 -2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,
239 -2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,
240 -2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,
241 -2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0,-2.0};
242 double * tgt_field = tgt_field_data;
243
244 if (!is_tgt)
245 for (size_t i = 0; i < grid_data.num_cells; ++i)
246 src_field[i] = (double)(grid_data.cell_ids[i] + 1);
247
248 yac_interpolation_execute(
249 interpolation, &src_fields, &tgt_field);
250
251 if (is_tgt) {
252
253 for (size_t i = 0; i < grid_data.num_vertices; ++i) {
254
255 double ref_tgt_value =
256 utest_compute_ref_tgt_value(
257 i, global_core_mask[1], global_field_mask[1],
258 global_core_mask[0], global_field_mask[0],
259 with_src_mask, with_tgt_mask,
260 weight_types[weight_types_idx], partial_coverage);
261
262 if (fabs(ref_tgt_value - tgt_field[i]) > 1e-3)
263 PUT_ERR("wrong interpolation result");
264 }
265 }
266 }
267
268 yac_interpolation_delete(interpolation);
269 } // reorder type
270 yac_interp_weights_delete(weights);
271 yac_interp_method_delete(method_stack);
272 } // partial coverage
273 } // weight type
274 yac_interp_grid_delete(interp_grid);
275 } // with tgt mask
276 } // with src mask
277 yac_dist_grid_pair_delete(grid_pair);
278 yac_basic_grid_delete(grids[1]);
279 yac_basic_grid_delete(grids[0]);
280 }
281
282 {// small test
283
284 // The global source grid has 5x4 cells:
285 //
286 //---------------
287 // setup
288 //---------------
289
290 int is_tgt = split_comm_size == 1;
291 double coordinates_x[2][4] =
292 {{-2,-1,1,2}, {-1.0,0.0,1.0}};
293 double coordinates_y[2][4] =
294 {{-2,-1,1,2}, {-1.0,0.0,1.0}};
295 double field_coordiantes_x[2][3*3] =
296 {{-1.5,0.0,1.5,
297 -1.5,0.0,1.5,
298 -1.5,0.0,1.5},
299 {-1.0,0.0,1.0,
300 -1.0,0.0,1.0,
301 -1.0,0.0,1.0}};
302 double field_coordiantes_y[2][3*3] =
303 {{-1.5,-1.5,-1.5,
304 0.0,0.0,0.0,
305 1.5,1.5,1.5},
306 {-1.0,-1.0,-1.0,
307 0.0,0.0,0.0,
308 1.0,1.0,1.0}};
309 int field_mask[3*3] = {0,1,1, 1,1,1, 1,1,0};
310 size_t const num_cells[2][2] = {{3,3}, {2,2}};
311 size_t local_start[2][2][2] = {{{0,0},{0,1}}, {{0,0}}};
312 size_t local_count[2][2][2] = {{{3,2},{3,2}}, {{2,2}}};
313 int with_halo = 0;
314 for (size_t i = 0; i <= num_cells[is_tgt][0]; ++i)
315 coordinates_x[is_tgt][i] *= YAC_RAD;
316 for (size_t i = 0; i <= num_cells[is_tgt][1]; ++i)
317 coordinates_y[is_tgt][i] *= YAC_RAD;
318
319 struct yac_basic_grid_data grid_data =
320 utest_generate_basic_grid_data_reg2d(
321 coordinates_x[is_tgt], coordinates_y[is_tgt], num_cells[is_tgt],
322 local_start[is_tgt][split_comm_rank],
323 local_count[is_tgt][split_comm_rank], with_halo);
324
325 struct yac_basic_grid * grids[2] =
326 {yac_basic_grid_new(grid_names[is_tgt], grid_data),
327 yac_basic_grid_empty_new(grid_names[is_tgt^1])};
328
329 double field_coords[3*3][3];
330 {
331 const yac_int * ids = (is_tgt)?grid_data.vertex_ids:grid_data.cell_ids;
332 size_t count = (is_tgt)?grid_data.num_vertices:grid_data.num_cells;
333 enum yac_location field_location[2] = {YAC_LOC_CELL, YAC_LOC_CORNER};
334 int src_field_mask[3*3];
335 for (size_t i = 0; i < count; ++i) {
336 LLtoXYZ_deg(
337 field_coordiantes_x[is_tgt][ids[i]],
338 field_coordiantes_y[is_tgt][ids[i]], field_coords[i]);
339 if (!is_tgt) src_field_mask[i] = field_mask[ids[i]];
340 }
341 yac_basic_grid_add_coordinates(
342 grids[0], field_location[is_tgt], field_coords, count);
343 if (!is_tgt)
344 yac_basic_grid_add_mask(
345 grids[0], field_location[is_tgt], src_field_mask, count, NULL);
346 }
347
348 struct yac_dist_grid_pair * grid_pair =
349 yac_dist_grid_pair_new(grids[0], grids[1], MPI_COMM_WORLD);
350
351 for (int with_src_mask = 0; with_src_mask < 2; ++with_src_mask) {
352
353 struct yac_interp_field src_fields[] =
354 {{.location = YAC_LOC_CELL,
355 .coordinates_idx = 0,
356 .masks_idx = with_src_mask?0:SIZE_MAX}};
357 size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
358 struct yac_interp_field tgt_field =
359 {.location = YAC_LOC_CORNER,
360 .coordinates_idx = 0,
361 .masks_idx = SIZE_MAX};
362
363 struct yac_interp_grid * interp_grid =
364 yac_interp_grid_new(grid_pair, grid_names[0], grid_names[1],
365 num_src_fields, src_fields, tgt_field);
366
367 for (int partial_coverage = 0; partial_coverage < 2; ++partial_coverage) {
368
369 struct interp_method * method_stack[] =
370 {yac_interp_method_ncc_new(
371 YAC_INTERP_NCC_DIST, partial_coverage),
372 yac_interp_method_fixed_new(-1.0), NULL};
373
374 struct yac_interp_weights * weights =
375 yac_interp_method_do_search(method_stack, interp_grid);
376
377 struct yac_interpolation * interpolation =
378 yac_interp_weights_get_interpolation(
379 weights, YAC_MAPPING_ON_SRC, 1,
380 YAC_FRAC_MASK_NO_VALUE, 1.0, 0.0, NULL, 1, 1);
381
382 // check generated interpolation
383 {
384 double src_field_data[3*3];
385 double * src_field = src_field_data;
386 double ** src_fields = &src_field;
387 double tgt_field_data[3*3] =
388 {-2.0,-2.0,-2.0,
389 -2.0,-2.0,-2.0,
390 -2.0,-2.0,-2.0};
391 double * tgt_field = tgt_field_data;
392
393 if (!is_tgt)
394 for (size_t i = 0; i < grid_data.num_cells; ++i)
395 src_field[i] = (double)(grid_data.cell_ids[i] + 1);
396
397 yac_interpolation_execute(
398 interpolation, &src_fields, &tgt_field);
399
400 if (is_tgt) {
401
402 yac_int tgt_to_src[2][2][9][4] =
403 {{{{0,1,3,4},{0,1,3,4},{1,2,4,5},
404 {0,1,3,4},{4},{1,2,4,5},
405 {3,4,6,7},{3,4,6,7},{4,5,7,8}},
406 {{0,1,3,4},{1,2,4,5},{1,2,4,5},
407 {3,4,6,7},{4},{4,5,7,8},
408 {3,4,6,7},{4,5,7,8},{4,5,7,8}}},
409 {{{1,3,4},{1,3,4},{1,2,4,5},
410 {1,3,4},{4},{1,2,4,5},
411 {3,4,6,7},{3,4,6,7},{4,5,7}},
412 {{1,3,4},{1,2,4,5},{1,2,4,5},
413 {3,4,6,7},{4},{4,5,7},
414 {3,4,6,7},{4,5,7},{4,5,7}}}};
415 int num_src_per_tgt[2][2][2][9] =
416 {{{{4,4,4, 4,1,4, 4,4,4},{4,4,4, 4,1,4, 4,4,4}},
417 {{0,0,4, 0,1,4, 4,4,0},{0,4,4, 4,1,0, 4,0,0}}},
418 {{{4,4,4, 4,1,4, 4,4,4},{4,4,4, 4,1,4, 4,4,4}},
419 {{3,3,4, 3,1,4, 4,4,3},{3,4,4, 4,1,3, 4,3,3}}}};
420
421 for (size_t i = 0; i < grid_data.num_vertices; ++i) {
422
423 yac_int tgt_id = grid_data.vertex_ids[i];
424 double ref_tgt_value[2] = {0.0, 0.0};
425
426 double tgt_coord[3];
427 LLtoXYZ_deg(
428 field_coordiantes_x[1][tgt_id], field_coordiantes_y[1][tgt_id],
429 tgt_coord);
430
431 for (int j = 0; j < 2; ++j) {
432
433 double inv_distance_sum = 0.0;
434 int num_src =
435 num_src_per_tgt[partial_coverage][with_src_mask][j][tgt_id];
436
437 if (num_src > 0) {
438 for (int k = 0; k < num_src; ++k) {
439 yac_int src_id = tgt_to_src[with_src_mask][j][tgt_id][k];
440 double src_coord[3];
441 LLtoXYZ_deg(
442 field_coordiantes_x[0][src_id],
443 field_coordiantes_y[0][src_id], src_coord);
444 double angle = get_vector_angle(src_coord, tgt_coord);
445 if (angle < yac_angle_tol) {
446 ref_tgt_value[j] = (double)(src_id + 1);
447 inv_distance_sum = 1.0;
448 break;
449 }
450 double inv_distance =
451 1.0 / get_vector_angle(src_coord, tgt_coord);
452 inv_distance_sum += inv_distance;
453 ref_tgt_value[j] += (double)(src_id + 1) * inv_distance;
454 }
455
456 ref_tgt_value[j] /= inv_distance_sum;
457 } else {
458 ref_tgt_value[j] = -1.0;
459 }
460 }
461
462 if ((fabs(ref_tgt_value[0] - tgt_field[i]) > 1e-3) &&
463 (fabs(ref_tgt_value[1] - tgt_field[i]) > 1e-3))
464 PUT_ERR("wrong interpolation result");
465 }
466 }
467 }
468
469 yac_interpolation_delete(interpolation);
470 yac_interp_weights_delete(weights);
471 yac_interp_method_delete(method_stack);
472 } // partial_coverage
473 yac_interp_grid_delete(interp_grid);
474 } // with_src_mask
475 yac_dist_grid_pair_delete(grid_pair);
476 yac_basic_grid_delete(grids[1]);
477 yac_basic_grid_delete(grids[0]);
478 }
479
480 MPI_Comm_free(&split_comm);
481
482 xt_finalize();
483
484 MPI_Finalize();
485
486 return TEST_EXIT_CODE;
487}
488
489static double utest_compute_ref_tgt_value(
490 size_t tgt_corner, int * tgt_global_core_mask, int * tgt_global_field_mask,
491 int * src_global_core_mask, int * src_global_field_mask, int with_src_mask,
492 int with_tgt_mask, enum yac_interp_ncc_weight_type weight_type,
493 int partial_coverage) {
494
495 // if target corner is masked
496 if (!tgt_global_core_mask[tgt_corner] ||
497 (with_tgt_mask && !tgt_global_field_mask[tgt_corner])) return -2.0;
498
499 // determine matching source cell
500 int source_grid_x_idx = ((int)tgt_corner % 10) / 2;
501 int source_grid_y_idx = ((int)tgt_corner / 10) / 2;
502 int source_cell = source_grid_x_idx + source_grid_y_idx * 5;
503
504 // check source cell core mask
505 if (!src_global_core_mask[source_cell]) return -1.0;
506
507 // determine closest source cell corner
508 int is_right_corner = (int)tgt_corner & 1;
509 int is_upper_corner = ((int)tgt_corner / 10) & 1;
510 int source_corner =
511 source_grid_x_idx + is_right_corner +
512 6 * (source_grid_y_idx + is_upper_corner);
513
514 // determine source cells surrounding the source vertex
515 int source_cells[4];
516 source_cells[0] = (5 * (source_corner / 6) + source_corner % 6) - 5 - 1;
517 source_cells[1] = (5 * (source_corner / 6) + source_corner % 6) - 5;
518 source_cells[2] = (5 * (source_corner / 6) + source_corner % 6) - 1;
519 source_cells[3] = (5 * (source_corner / 6) + source_corner % 6);
520
521 // special handling for grid edge corners
522 if (source_corner <= 5) {
523 source_cells[0] = -1;
524 source_cells[1] = -1;
525 }
526 if (source_corner >= 24) {
527 source_cells[2] = -1;
528 source_cells[3] = -1;
529 }
530 if ((source_corner % 6) == 0) {
531 source_cells[0] = -1;
532 source_cells[2] = -1;
533 }
534 if (((source_corner + 1) % 6) == 0) {
535 source_cells[3] = -1;
536 source_cells[1] = -1;
537 }
538
539 // apply core mask to source cells
540 for (int i = 0; i < 4; ++i)
541 if ((source_cells[i] >= 0) &&
542 (!src_global_core_mask[source_cells[i]]))
543 source_cells[i] = -1;
544
545 // apply field mask to source cells
546 if (with_src_mask) {
547 for (int i = 0; i < 4; ++i) {
548 if ((source_cells[i] >= 0) &&
549 (!src_global_field_mask[source_cells[i]])) {
550 if (!partial_coverage) return -1.0;
551 source_cells[i] = -1;
552 }
553 }
554 }
555
556 double const ref_distances[2][2] =
557 {{sqrt(0.25*0.25 + 0.25*0.25),
558 sqrt(0.75*0.75 + 0.25*0.25)},
559 {sqrt(0.25*0.25 + 0.75*0.75),
560 sqrt(0.75*0.75 + 0.75*0.75)}};
561
562 double distances[4];
563 int num_src_cells = 0;
564 for (int i = 0, k = 0; i < 2; ++i) {
565 for (int j = 0; j < 2; ++j, ++k) {
566 if (source_cells[k] >= 0) {
567 distances[num_src_cells] =
568 ref_distances[i^is_upper_corner^1][j^is_right_corner^1];
569 source_cells[num_src_cells] = source_cells[k];
570 ++num_src_cells;
571 }
572 }
573 }
574
575 if (num_src_cells == 0) return -1.0;
576
577 double tgt_value = 0.0;
578 switch(weight_type) {
579 default:
580 case (YAC_INTERP_NCC_AVG): {
581 double weight = 1.0 / (double)num_src_cells;
582 for (int i = 0; i < num_src_cells; ++i)
583 tgt_value += (double)(source_cells[i] + 1) * weight;
584 break;
585 }
586 case (YAC_INTERP_NCC_DIST): {
587 double inv_distance_sum = 0.0;
588 for (int i = 0; i < num_src_cells; ++i)
589 inv_distance_sum += 1.0 / distances[i];
590 for (int i = 0; i < num_src_cells; ++i)
591 tgt_value +=
592 (double)(source_cells[i] + 1) / (distances[i] * inv_distance_sum);
593 break;
594 }
595 };
596 return tgt_value;
597}
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
size_t yac_basic_grid_add_coordinates(struct yac_basic_grid *grid, enum yac_location location, yac_coordinate_pointer coordinates, size_t count)
Definition basic_grid.c:222
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_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
yac_angle_tol
#define yac_angle_tol
Definition geometry.h:26
get_vector_angle
static double get_vector_angle(double const a[3], double const b[3])
Definition geometry.h:472
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_fixed_new
struct interp_method * yac_interp_method_fixed_new(double value)
Definition interp_method_fixed.c:63
interp_method_fixed.h
yac_interp_method_ncc_new
struct interp_method * yac_interp_method_ncc_new(enum yac_interp_ncc_weight_type weight_type, int partial_coverage)
Definition interp_method_ncc.c:432
interp_method_ncc.h
yac_interp_ncc_weight_type
yac_interp_ncc_weight_type
Definition interp_method_ncc.h:12
YAC_INTERP_NCC_DIST
@ YAC_INTERP_NCC_DIST
distance weighted average of n source points
Definition interp_method_ncc.h:14
YAC_INTERP_NCC_AVG
@ YAC_INTERP_NCC_AVG
average of n source points
Definition interp_method_ncc.h:13
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_location
yac_location
Definition location.h:12
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_plot_grids.ids
ids
Definition yac_plot_grids.py:136
grid_data
Definition test_dummy_coupling3_c.c:36
interp_method
Definition interp_method_internal.h:36
yac_basic_grid_data
Definition basic_grid_data.h:13
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
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
grid_names
static char const * grid_names[2]
Definition test_interp_method_ncc_parallel.c:28
num_reorder_types
size_t num_reorder_types
Definition test_interp_method_ncc_parallel.c:27
reorder_types
enum yac_interp_weights_reorder_type reorder_types[]
Definition test_interp_method_ncc_parallel.c:25
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
weight_type
weight_type
Definition test_interp_weights_parallel.c:32
coordinates_x
double coordinates_x[]
Definition test_interpolation_parallel3.c:105
num_cells
size_t num_cells[2]
Definition test_interpolation_parallel3.c:104
coordinates_y
double coordinates_y[]
Definition test_interpolation_parallel3.c:106
with_halo
int with_halo
Definition test_interpolation_parallel3.c:110
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