18#define AREA_TOL_FACTOR (1e-6)
22 size_t * tgt_points,
size_t count,
26 size_t * tgt_points,
size_t count,
79 int max_num_vertices_per_cell = 0;
82 max_num_vertices_per_cell)
83 max_num_vertices_per_cell =
85 return max_num_vertices_per_cell;
97 *src_grid_cells =
xmalloc(max_num_src_per_tgt *
sizeof(**src_grid_cells));
99 double (*coordinates_xyz_buffer)[3];
103 int max_num_vertices_per_cell =
108 xmalloc((max_num_src_per_tgt + 1) * (
size_t)max_num_vertices_per_cell *
109 sizeof(*edge_type_buffer));
110 coordinates_xyz_buffer =
111 xmalloc((max_num_src_per_tgt + 1) * (
size_t)max_num_vertices_per_cell *
112 sizeof(*coordinates_xyz_buffer));
115 tgt_grid_cell->
edge_type = edge_type_buffer;
116 tgt_grid_cell->
array_size = max_num_vertices_per_cell;
117 for (
size_t i = 0; i < max_num_src_per_tgt; ++i) {
118 (*src_grid_cells)[i].coordinates_xyz =
119 coordinates_xyz_buffer + (i + 1) * max_num_vertices_per_cell;
120 (*src_grid_cells)[i].edge_type =
121 edge_type_buffer + (i + 1) * max_num_vertices_per_cell;
122 (*src_grid_cells)[i].array_size = max_num_vertices_per_cell;
136 int max_num_vertices_per_cell =
141 xmalloc(2 * (
size_t)max_num_vertices_per_cell *
142 sizeof(*edge_type_buffer));
144 xmalloc(2 * (
size_t)max_num_vertices_per_cell *
145 sizeof(*coordinates_xyz_buffer));
148 tgt_grid_cell->
edge_type = edge_type_buffer;
149 tgt_grid_cell->
array_size = max_num_vertices_per_cell;
152 coordinates_xyz_buffer + max_num_vertices_per_cell;
154 edge_type_buffer + max_num_vertices_per_cell;
155 src_grid_cell->
array_size = max_num_vertices_per_cell;
161 size_t * src_cells,
struct yac_grid_cell tgt_grid_cell_buffer,
163 double * weights,
size_t * num_weights,
int partial_coverage,
165 int enforced_conserv) {
168 tgt_basic_grid_data, tgt_cell, &tgt_grid_cell_buffer);
169 for (
size_t i = 0; i < src_count; ++i)
171 src_basic_grid_data, src_cells[i], src_grid_cell_buffer + i);
173 double * area = weights;
175 src_count, src_grid_cell_buffer, tgt_grid_cell_buffer, area);
177 size_t num_valid_weights = 0;
178 for (
size_t i = 0; i < src_count; ++i) {
181 if (i != num_valid_weights) {
182 area[num_valid_weights] = area[i];
183 src_cells[num_valid_weights] = src_cells[i];
188 *num_weights = num_valid_weights;
189 if (num_valid_weights == 0)
return 0;
197 "ERROR(compute_weights_order_first_conserv_no_partial): "
198 "invalid normalisation option in conservative remapping")
199 switch(normalisation) {
201 norm_factor = 1.0 / tgt_cell_area;
205 double fracarea = 0.0;
206 for (
size_t i = 0; i < num_valid_weights; ++i) fracarea += area[i];
207 norm_factor = 1.0 / fracarea;
212 if (partial_coverage) {
213 for (
size_t i = 0; i < num_valid_weights; ++i) weights[i] *= norm_factor;
216 double tgt_cell_area_diff = tgt_cell_area;
218 for (
size_t i = 0; i < num_valid_weights; ++i) {
219 double curr_area = area[i];
220 tgt_cell_area_diff -= curr_area;
221 weights[i] = curr_area * norm_factor;
223 int successful = fabs(tgt_cell_area_diff) <= area_tol;
224 if (successful && enforced_conserv)
232 size_t * tgt_points,
size_t count,
240 "ERROR(do_search_conserv): invalid number of source fields")
244 "ERROR(do_search_conserv): unsupported source field location type")
248 "ERROR(do_search_conserv): unsupported target field location type")
251 size_t * src_cells = NULL;
252 size_t * num_src_per_tgt =
xmalloc(count *
sizeof(*num_src_per_tgt));
256 interp_grid, tgt_points, count, &src_cells, num_src_per_tgt);
265 size_t total_num_weights = 0;
266 size_t max_num_src_per_tgt = 0;
267 for (
size_t i = 0; i <
count; ++i) {
268 size_t curr_num_src_per_tgt = num_src_per_tgt[i];
269 if (curr_num_src_per_tgt > max_num_src_per_tgt)
270 max_num_src_per_tgt = curr_num_src_per_tgt;
271 total_num_weights += num_src_per_tgt[i];
277 yac_int * temp_src_global_ids =
278 xmalloc(max_num_src_per_tgt *
sizeof(*temp_src_global_ids));
280 for (
size_t i = 0, offset = 0; i <
count; ++i) {
282 size_t curr_num_src_per_tgt = num_src_per_tgt[i];
283 size_t * curr_src_cells = src_cells + offset;
284 offset += curr_num_src_per_tgt;
286 for (
size_t j = 0; j < curr_num_src_per_tgt; ++j)
287 temp_src_global_ids[j] =
291 temp_src_global_ids, curr_num_src_per_tgt, curr_src_cells);
294 free(temp_src_global_ids);
297 double * w =
xmalloc(total_num_weights *
sizeof(*w));
298 size_t result_count = 0;
299 size_t * failed_tgt =
xmalloc(
count *
sizeof(*failed_tgt));
300 total_num_weights = 0;
310 interp_grid, max_num_src_per_tgt, &tgt_grid_cell, &src_grid_cells);
313 for (
size_t i = 0, offset = 0, result_offset = 0; i < count; ++i) {
315 size_t curr_src_count = num_src_per_tgt[i];
316 size_t curr_tgt_point = tgt_points[i];
321 tgt_basic_grid_data, curr_tgt_point,
322 src_basic_grid_data, curr_src_count, src_cells + offset,
323 tgt_grid_cell, src_grid_cells, w + result_offset, &num_weights,
324 partial_coverage, normalisation, enforced_conserv)) {
326 if (offset != result_offset) {
329 src_cells + result_offset, src_cells + offset,
330 num_weights *
sizeof(*src_cells));
332 tgt_points[result_count] = curr_tgt_point;
333 num_src_per_tgt[result_count] = num_weights;
335 result_offset += num_weights;
336 total_num_weights += num_weights;
338 failed_tgt[i - result_count] = curr_tgt_point;
341 offset += curr_src_count;
346 free(src_grid_cells);
348 if (result_count != count)
349 memcpy(tgt_points + result_count, failed_tgt,
350 (count - result_count) *
sizeof(*tgt_points));
356 interp_grid, tgt_points, result_count),
357 .count = result_count};
360 interp_grid, 0, src_cells, total_num_weights);
364 weights, &tgts, num_src_per_tgt, srcs, w);
369 free(num_src_per_tgt);
420 for (
size_t k = 0; k < 3; ++k)
421 for (
size_t l = 0; l < 3; ++l)
422 M[k][l] = - src_cell_centroid[k] * src_cell_centroid[l];
423 for (
size_t k = 0; k < 3; ++k)
429 for (
size_t i = 0; i < N; ++i) {
432 for (
size_t j = 0; j < 3; ++j) buffer[i].
weight[j] = 0.0;
435 for (
size_t n = 0; n < N; ++n)
436 for (
size_t i = 0; i < 3; ++i)
437 for (
size_t j = 0; j < 3; ++j)
440 memcpy(G_i->
data, buffer, N *
sizeof(*buffer));
444 void const * a,
void const * b) {
453 double abs_weight_a = fabs(weight_a->
weight);
454 double abs_weight_b = fabs(weight_b->
weight);
455 ret = (abs_weight_a > abs_weight_b) - (abs_weight_a < abs_weight_b);
462 void const * a,
void const * b) {
485 for (
size_t i = 1; i < n_; ++i, ++curr_weight_data) {
493 *prev_weight_data = *curr_weight_data;
501 for (
size_t i = 0; i < n_; ++i) {
503 if (weights[i].
weight == 0.0)
continue;
504 if (i != new_n) weights[new_n] = weights[i];
520 size_t N = src_cell_gradient->
n;
525 double * overlap_barycenter = super_cell->
barycenter;
527 for (
size_t n = 1; n <= N; ++n) {
529 (gradient_weights[n].
weight[0] * overlap_barycenter[0] +
530 gradient_weights[n].
weight[1] * overlap_barycenter[1] +
531 gradient_weights[n].
weight[2] * overlap_barycenter[2]) *
543 double barycenter[3]) {
546 size_t const * vertices =
553 for (
size_t i = 0; i < num_vertices; ++i) {
554 double const * curr_vertex_coordinate =
556 barycenter[0] += curr_vertex_coordinate[0];
557 barycenter[1] += curr_vertex_coordinate[1];
558 barycenter[2] += curr_vertex_coordinate[2];
564 struct yac_interp_grid * interp_grid,
size_t * tgt_points,
size_t count,
566 int * interp_fail_flag,
size_t ** src_cells,
size_t * num_src_cells,
568 int partial_coverage) {
573 "ERROR(compute_super_cells): "
574 "invalid normalisation option in conservative remapping")
576 size_t * num_src_per_tgt =
xmalloc(count *
sizeof(*num_src_per_tgt));
580 interp_grid, tgt_points, count, src_cells, num_src_per_tgt);
583 size_t total_num_overlaps = 0;
584 for (
size_t i = 0; i < count; ++i) total_num_overlaps += num_src_per_tgt[i];
586 *num_src_cells = total_num_overlaps;
589 size_t * num_tgt_per_src =
590 xrealloc(num_src_per_tgt, *num_src_cells *
sizeof(*num_tgt_per_src));
595 size_t * tgt_cells = NULL;
597 interp_grid, *src_cells, *num_src_cells, &tgt_cells, num_tgt_per_src);
599 total_num_overlaps = 0;
600 for (
size_t i = 0; i < *num_src_cells; ++i)
601 total_num_overlaps += num_tgt_per_src[i];
604 xmalloc(total_num_overlaps *
sizeof(*super_cells));
611 for (
size_t i = 0, j = 0; i < *num_src_cells; ++i) {
613 size_t curr_num_overlaps = num_tgt_per_src[i];
614 size_t curr_src_cell = (*src_cells)[i];
618 for (
size_t k = 0; k < curr_num_overlaps; ++k, ++j) {
620 size_t curr_tgt_cell = tgt_cells[j];
630 free(num_tgt_per_src);
634 qsort(super_cells, total_num_overlaps,
sizeof(*super_cells),
646 size_t new_num_super_cells = 0;
648 for (
size_t i = 0, j = 0; i < total_num_overlaps;) {
651 size_t curr_tgt_cell = super_cells[i].
tgt.
local_id;
653 tgt_basic_grid_data, curr_tgt_cell, &tgt_grid_cell);
654 double curr_tgt_cell_coverage = 0.0;
657 for (;(i < total_num_overlaps) &&
658 (super_cells[i].tgt.
local_id == curr_tgt_cell); ++i) {
662 src_basic_grid_data, super_cells[i].
src.local_id, &src_grid_cell);
665 double super_cell_area;
666 double barycenter[3];
668 1, &src_grid_cell, tgt_grid_cell, &super_cell_area, &barycenter);
671 if (super_cell_area > 0.0) {
673 super_cells[new_num_super_cells].
src = super_cells[i].
src;
674 super_cells[new_num_super_cells].
tgt = super_cells[i].
tgt;
675 super_cells[new_num_super_cells].
area = super_cell_area;
676 memcpy(super_cells[new_num_super_cells].barycenter, barycenter,
679 ++new_num_super_cells;
681 curr_tgt_cell_coverage += super_cell_area;
688 if (new_num_super_cells != j) {
694 "ERROR(compute_super_cells): invalid normalisation")
695 switch (normalisation) {
698 norm_factor = 1.0 / curr_tgt_cell_area;
701 norm_factor = 1.0 / curr_tgt_cell_coverage;
705 for (; j < new_num_super_cells; ++j)
706 super_cells[j].norm_area = super_cells[j].area * norm_factor;
711 while ((tgt_idx < count) && (tgt_points[tgt_idx] < curr_tgt_cell))
712 interp_fail_flag[tgt_idx++] = 1;
714 if ((tgt_idx < count) && (tgt_points[tgt_idx] == curr_tgt_cell)) {
718 if (partial_coverage) {
719 interp_fail_flag[tgt_idx] = curr_tgt_cell_coverage < area_tol;
721 interp_fail_flag[tgt_idx] =
722 fabs(curr_tgt_cell_area - curr_tgt_cell_coverage) > area_tol;
728 for (; tgt_idx < count; ++tgt_idx) interp_fail_flag[tgt_idx] = 1;
729 *num_super_cells = new_num_super_cells;
731 xrealloc(super_cells, new_num_super_cells *
sizeof(*super_cells));
738 size_t * src_cells,
int * skip_src_cell,
size_t num_src_cells,
742 xmalloc(num_src_cells *
sizeof(*src_cell_centroids));
746 qsort(super_cells, num_super_cells,
sizeof(*super_cells),
762 for (
size_t i = 0, offset = 0; i < num_src_cells; ++i) {
764 if (skip_src_cell[i])
continue;
766 size_t curr_src_cell = src_cells[i];
768 src_basic_grid_data, curr_src_cell, &src_grid_cell);
772 size_t curr_num_super_cells = offset;
773 while ((offset < num_super_cells) &&
774 (super_cells[offset].
src.local_id == curr_src_cell)) ++offset;
775 curr_num_super_cells = offset - curr_num_super_cells;
777 double src_cell_area_diff = src_cell_area;
778 double src_cell_centroid[3] = {0.0, 0.0, 0.0};
780 for (
size_t j = 0; j < curr_num_super_cells; ++j) {
782 double super_cell_area = curr_super_cells[j].
area;
783 double * super_cell_barycenter = curr_super_cells[j].
barycenter;
784 src_cell_centroid[0] += super_cell_area * super_cell_barycenter[0];
785 src_cell_centroid[1] += super_cell_area * super_cell_barycenter[1];
786 src_cell_centroid[2] += super_cell_area * super_cell_barycenter[2];
787 src_cell_area_diff -= super_cell_area;
791 src_cell_centroids[i][0] = src_cell_centroid[0];
792 src_cell_centroids[i][1] = src_cell_centroid[1];
793 src_cell_centroids[i][2] = src_cell_centroid[2];
799 return src_cell_centroids;
805 size_t num_src_cells,
size_t * src_cell_neighbours) {
808 xmalloc(num_src_cells *
sizeof(*src_cell_gradients));
813 size_t total_num_gradient_weights = 0;
814 size_t max_num_neigh_per_src = 0;
815 for (
size_t i = 0; i < num_src_cells; ++i) {
816 src_cell_gradients[i].
data = NULL;
817 src_cell_gradients[i].
n = 0;
818 size_t curr_num_neigh =
820 total_num_gradient_weights += curr_num_neigh + 1;
821 if (max_num_neigh_per_src < curr_num_neigh)
822 max_num_neigh_per_src = curr_num_neigh;
825 (total_num_gradient_weights > 0)?
826 (
xmalloc(total_num_gradient_weights *
827 sizeof(*weight_vector_data_buffer))):NULL;
828 for (
size_t i = 0; i < total_num_gradient_weights; ++i)
829 for (
size_t j = 0; j < 3; ++j)
830 weight_vector_data_buffer[i].
weight[j] = 0.0;
833 xmalloc((max_num_neigh_per_src + 1) *
sizeof(*orth_buffer));
857 for (
size_t i = 0, offset = 0, weight_vector_data_buffer_offset = 0;
858 i < num_src_cells; ++i) {
860 size_t curr_num_neigh =
862 size_t * curr_neighs = src_cell_neighbours + offset;
863 offset += curr_num_neigh;
865 G_i->
data = weight_vector_data_buffer + weight_vector_data_buffer_offset;
867 if (skip_src_cell[i])
continue;
869 weight_vector_data_buffer_offset += curr_num_neigh + 1;
870 G_i->
n = curr_num_neigh + 1;
874 for (
size_t j = 0; j < curr_num_neigh; ++j) {
875 size_t curr_neigh = curr_neighs[j];
877 if (curr_neigh != SIZE_MAX) {
896 .num_corners = 3, .array_size = 0};
906 for (
size_t j = 0; j < curr_num_neigh; ++j) {
908 size_t neigh_idx[2] = {j + 1, (j+1)%curr_num_neigh+1};
909 size_t neigh_local_ids[2] =
913 if (neigh_local_ids[0] == neigh_local_ids[1])
continue;
915 double neigh_cell_barycenters[2][3];
917 src_basic_grid_data, neigh_local_ids[0], neigh_cell_barycenters[0]);
919 src_basic_grid_data, neigh_local_ids[1], neigh_cell_barycenters[1]);
921 centroid_triangle.
coordinates_xyz[1][0] = neigh_cell_barycenters[0][0];
922 centroid_triangle.
coordinates_xyz[1][1] = neigh_cell_barycenters[0][1];
923 centroid_triangle.
coordinates_xyz[1][2] = neigh_cell_barycenters[0][2];
924 centroid_triangle.
coordinates_xyz[2][0] = neigh_cell_barycenters[1][0];
925 centroid_triangle.
coordinates_xyz[2][1] = neigh_cell_barycenters[1][1];
926 centroid_triangle.
coordinates_xyz[2][2] = neigh_cell_barycenters[1][2];
933 neigh_cell_barycenters[0], neigh_cell_barycenters[1], C_j_x_C_k);
935 double curr_edge_direction =
948 for (
size_t l = 0; l < 3; ++l) C_j_x_C_k[l] *= 0.5;
951 G_i->
data[neigh_idx[0]].
weight[0] += C_j_x_C_k[0];
952 G_i->
data[neigh_idx[0]].
weight[1] += C_j_x_C_k[1];
953 G_i->
data[neigh_idx[0]].
weight[2] += C_j_x_C_k[2];
956 G_i->
data[neigh_idx[1]].
weight[0] += C_j_x_C_k[0];
957 G_i->
data[neigh_idx[1]].
weight[1] += C_j_x_C_k[1];
958 G_i->
data[neigh_idx[1]].
weight[2] += C_j_x_C_k[2];
963 for (
size_t j = 0; j <= curr_num_neigh; ++j)
964 for (
size_t k = 0; k < 3; ++k)
967 double inv_A_C_K = (A_C_k >
YAC_AREA_TOL)?(1.0/A_C_k):0.0;
969 for (
size_t k = 0; k <= curr_num_neigh; ++k)
970 for (
size_t l = 0; l < 3; ++l)
974 src_cell_centroids[i], G_i, orth_buffer);
978 return src_cell_gradients;
983 int * interp_fail_flag,
size_t num_tgt_cells,
985 size_t ** src_per_tgt,
double ** weights,
size_t * num_src_per_tgt) {
987 size_t num_interpolated_tgt = 0;
991 qsort(super_cells, num_super_cells,
sizeof(*super_cells),
994 size_t max_num_weights_per_tgt = 0;
995 size_t max_num_total_weights = 0;
998 for (
size_t i = 0, j = 0; i < num_tgt_cells; ++i) {
1000 if (interp_fail_flag[i])
continue;
1002 size_t curr_tgt_cell = tgt_cells[i];
1003 size_t curr_num_weights = 0;
1006 while ((j < num_super_cells) &&
1007 (super_cells[j].tgt.
local_id < curr_tgt_cell)) ++j;
1010 while ((j < num_super_cells) &&
1011 (super_cells[j].tgt.
local_id == curr_tgt_cell))
1014 max_num_total_weights += curr_num_weights;
1015 if (max_num_weights_per_tgt < curr_num_weights)
1016 max_num_weights_per_tgt = curr_num_weights;
1017 num_interpolated_tgt++;
1021 xmalloc(max_num_weights_per_tgt *
sizeof(*weight_buffer));
1022 *weights =
xmalloc(max_num_total_weights *
sizeof(**weights));
1023 *src_per_tgt =
xmalloc(max_num_total_weights *
sizeof(**src_per_tgt));
1044 for (
size_t i = 0, j = 0; i < num_interpolated_tgt; ++i) {
1046 size_t curr_tgt_cell = tgt_cells[i];
1049 while ((j < num_super_cells) &&
1050 (super_cells[j].tgt.
local_id != curr_tgt_cell)) ++j;
1052 size_t weight_buffer_offset = 0;
1055 while ((j < num_super_cells) &&
1056 (super_cells[j].tgt.
local_id == curr_tgt_cell)) {
1058 size_t num_weights =
1060 super_cells + j, weight_buffer + weight_buffer_offset);
1061 weight_buffer_offset += num_weights;
1068 num_src_per_tgt[i] = weight_buffer_offset;
1070 for (
size_t k = 0; k < weight_buffer_offset; ++k, ++w_idx) {
1071 (*weights)[w_idx] = weight_buffer[k].
weight;
1072 (*src_per_tgt)[w_idx] = weight_buffer[k].
local_id;
1075 free(weight_buffer);
1077 return num_interpolated_tgt;
1082 size_t * tgt_points,
size_t count,
1090 "ERROR(do_search_conserv): invalid number of source fields")
1094 "ERROR(do_search_conserv): unsupported source field location type")
1098 "ERROR(do_search_conserv): unsupported target field location type")
1103 int * interp_fail_flag =
xmalloc(count *
sizeof(*interp_fail_flag));
1105 size_t num_super_cells = 0;
1106 size_t * src_cells = NULL;
1107 size_t num_src_cells = 0;
1111 interp_grid, tgt_points, count, &super_cells, &num_super_cells,
1112 interp_fail_flag, &src_cells, &num_src_cells,
1115 size_t total_num_src_cell_neighbours = 0;
1118 for (
size_t i = 0; i < num_src_cells; ++i)
1119 total_num_src_cell_neighbours +=
1121 size_t * src_cell_neighbours =
1122 xmalloc(total_num_src_cell_neighbours *
sizeof(*src_cell_neighbours));
1126 interp_grid, src_cells, num_src_cells, src_cell_neighbours);
1130 int * skip_src_cell =
xmalloc(num_src_cells *
sizeof(*skip_src_cell));
1134 if (src_cell_mask != NULL) {
1135 for (
size_t i = 0, offset = 0; i < num_src_cells; ++i) {
1136 size_t curr_src_cell = src_cells[i];
1137 if ((skip_src_cell[i] = !src_cell_mask[curr_src_cell]))
continue;
1138 size_t * curr_neighbours = src_cell_neighbours + offset;
1139 size_t curr_num_neigh =
1141 offset += curr_num_neigh;
1142 for (
size_t j = 0; j < curr_num_neigh; ++j)
1143 if ((curr_neighbours[j] != SIZE_MAX) &&
1144 (!src_cell_mask[curr_neighbours[j]]))
1145 curr_neighbours[j] = SIZE_MAX;
1148 memset(skip_src_cell, 0, num_src_cells *
sizeof(*skip_src_cell));
1154 num_src_cells, super_cells, num_super_cells);
1159 interp_grid, src_cells, src_cell_centroids,
1160 skip_src_cell, num_src_cells, src_cell_neighbours);
1161 free(src_cell_neighbours);
1162 free(src_cell_centroids);
1165 qsort(super_cells, num_super_cells,
sizeof(*super_cells),
1167 for (
size_t i = 0, j = 0; i < num_src_cells; ++i) {
1168 size_t curr_src_cell = src_cells[i];
1173 (j >= num_super_cells) ||
1174 (super_cells[j].
src.local_id >= curr_src_cell),
1175 "ERROR(do_search_conserv_2nd_order): internal error");
1176 while ((j < num_super_cells) &&
1177 (super_cells[j].
src.local_id == curr_src_cell)) {
1181 free(skip_src_cell);
1184 size_t * src_per_tgt = NULL;
1186 size_t * num_src_per_tgt =
xmalloc(count *
sizeof(*num_src_per_tgt));
1189 size_t num_interpolated_tgt =
1191 interp_grid, tgt_points, interp_fail_flag, count,
1192 super_cells, num_super_cells, &src_per_tgt, &w, num_src_per_tgt);
1193 if (num_src_cells > 0) free(src_cell_gradients->
data);
1194 free(src_cell_gradients);
1196 free(interp_fail_flag);
1198 size_t total_num_weights = 0;
1199 for (
size_t i = 0; i < num_interpolated_tgt; ++i)
1200 total_num_weights += num_src_per_tgt[i];
1205 interp_grid, tgt_points, num_interpolated_tgt),
1206 .count = num_interpolated_tgt};
1209 interp_grid, 0, src_per_tgt, total_num_weights);
1213 weights, &tgts, num_src_per_tgt, srcs, w);
1218 free(num_src_per_tgt);
1221 return num_interpolated_tgt;
1225 int order,
int enforced_conserv,
int partial_coverage,
1232 "ERROR(yac_interp_method_conserv_new): interp_method_conserv only "
1233 "supports enforced_conserv with first order conservative remapping")
1235 (order == 1) || (order == 2),
1236 "ERROR(yac_interp_method_conserv_new): invalid order")
double yac_huiliers_area(struct yac_grid_cell cell)
Area calculation on a unit sphere taken from ESMF based on L'Huilier's Theorem.
Structs and interfaces for area calculations.
static int edge_direction(double *a, double *b)
void yac_correct_weights(size_t nSourceCells, double *weight)
correct interpolation weights
void yac_compute_overlap_info(size_t N, struct yac_grid_cell *source_cell, struct yac_grid_cell target_cell, double *overlap_areas, double(*overlap_barycenters)[3])
calculates partial areas for all overlapping parts of the source cells with arbitrary target cells,...
void yac_compute_overlap_areas(size_t N, struct yac_grid_cell *source_cell, struct yac_grid_cell target_cell, double *partial_areas)
calculates partial areas for all overlapping parts of the source cells with arbitrary target cells,...
void yac_const_basic_grid_data_get_grid_cell(struct yac_const_basic_grid_data *grid_data, size_t cell_idx, struct yac_grid_cell *buffer_cell)
int const *const const_int_pointer
static void crossproduct_kahan(double const a[], double const b[], double cross[])
static void normalise_vector(double v[])
@ YAC_GREAT_CIRCLE_EDGE
great circle
void yac_interp_grid_do_cell_search_tgt(struct yac_interp_grid *interp_grid, size_t *src_cells, size_t count, size_t **tgt_cells, size_t *num_tgt_per_src)
enum yac_location yac_interp_grid_get_tgt_field_location(struct yac_interp_grid *interp_grid)
size_t yac_interp_grid_get_num_src_fields(struct yac_interp_grid *interp_grid)
const_int_pointer yac_interp_grid_get_src_field_mask(struct yac_interp_grid *interp_grid, size_t src_field_idx)
void yac_interp_grid_get_src_cell_neighbours(struct yac_interp_grid *interp_grid, size_t *src_cells, size_t count, size_t *neighbours)
struct remote_point * yac_interp_grid_get_tgt_remote_points(struct yac_interp_grid *interp_grid, size_t *tgt_points, size_t count)
enum yac_location yac_interp_grid_get_src_field_location(struct yac_interp_grid *interp_grid, size_t src_field_idx)
struct remote_point * yac_interp_grid_get_src_remote_points(struct yac_interp_grid *interp_grid, size_t src_field_idx, size_t *src_points, size_t count)
struct yac_const_basic_grid_data * yac_interp_grid_get_basic_grid_data_tgt(struct yac_interp_grid *interp_grid)
void yac_interp_grid_do_cell_search_src(struct yac_interp_grid *interp_grid, size_t *tgt_cells, size_t count, size_t **src_cells, size_t *num_src_per_tgt)
struct yac_const_basic_grid_data * yac_interp_grid_get_basic_grid_data_src(struct yac_interp_grid *interp_grid)
static void delete_conserv(struct interp_method *method)
static void get_cell_buffers_(struct yac_interp_grid *interp_grid, struct yac_grid_cell *tgt_grid_cell, struct yac_grid_cell *src_grid_cell)
static size_t compute_2nd_order_weights(struct yac_interp_grid *interp_grid, size_t *tgt_cells, int *interp_fail_flag, size_t num_tgt_cells, struct supermesh_cell *super_cells, size_t num_super_cells, size_t **src_per_tgt, double **weights, size_t *num_src_per_tgt)
static yac_coordinate_pointer compute_src_cell_centroids(struct yac_interp_grid *interp_grid, size_t *src_cells, int *skip_src_cell, size_t num_src_cells, struct supermesh_cell *super_cells, size_t num_super_cells)
static void get_cell_buffers(struct yac_interp_grid *interp_grid, size_t max_num_src_per_tgt, struct yac_grid_cell *tgt_grid_cell, struct yac_grid_cell **src_grid_cells)
static void compute_cell_barycenter(struct yac_const_basic_grid_data *grid_data, size_t cell_idx, double barycenter[3])
static size_t do_search_conserv_1st_order(struct interp_method *method, struct yac_interp_grid *interp_grid, size_t *tgt_points, size_t count, struct yac_interp_weights *weights)
static void compact_weight_vector_data(struct weight_vector_data *weights, size_t *n)
static struct interp_method_vtable interp_method_conserv_2nd_order_vtable
struct interp_method * yac_interp_method_conserv_new(int order, int enforced_conserv, int partial_coverage, enum yac_interp_method_conserv_normalisation normalisation)
static int get_max_num_vertices_per_cell(struct yac_const_basic_grid_data *basic_grid_data)
static void compute_super_cells(struct yac_interp_grid *interp_grid, size_t *tgt_points, size_t count, struct supermesh_cell **super_cells_, size_t *num_super_cells, int *interp_fail_flag, size_t **src_cells, size_t *num_src_cells, enum yac_interp_method_conserv_normalisation normalisation, int partial_coverage)
static int compare_supermesh_cell_tgt_local_ids(const void *a, const void *b)
static void orthogonalise_weight_vector(double *src_cell_centroid, struct weight_vector_3d *G_i, struct weight_vector_data_3d *buffer)
static struct interp_method_vtable interp_method_conserv_1st_order_vtable
static int compare_supermesh_cell_src_local_ids(const void *a, const void *b)
static size_t compute_2nd_order_tgt_cell_weights(struct supermesh_cell *super_cell, struct weight_vector_data *weights)
static int compare_weight_vector_data_weight(void const *a, void const *b)
static int compare_weight_vector_data(void const *a, void const *b)
static size_t do_search_conserv_2nd_order(struct interp_method *method, struct yac_interp_grid *interp_grid, size_t *tgt_points, size_t count, struct yac_interp_weights *weights)
static struct weight_vector_3d * compute_src_cell_gradients(struct yac_interp_grid *interp_grid, size_t *src_cells, yac_coordinate_pointer src_cell_centroids, int *skip_src_cell, size_t num_src_cells, size_t *src_cell_neighbours)
static int compute_1st_order_weights(struct yac_const_basic_grid_data *tgt_basic_grid_data, size_t tgt_cell, struct yac_const_basic_grid_data *src_basic_grid_data, size_t src_count, size_t *src_cells, struct yac_grid_cell tgt_grid_cell_buffer, struct yac_grid_cell *src_grid_cell_buffer, double *weights, size_t *num_weights, int partial_coverage, enum yac_interp_method_conserv_normalisation normalisation, int enforced_conserv)
yac_interp_method_conserv_normalisation
@ YAC_INTERP_CONSERV_DESTAREA
@ YAC_INTERP_CONSERV_FRACAREA
void yac_interp_weights_add_wsum(struct yac_interp_weights *weights, struct remote_points *tgts, size_t *num_src_per_tgt, struct remote_point *srcs, double *w)
#define xrealloc(ptr, size)
enum yac_interp_method_conserv_normalisation normalisation
struct interp_method_vtable * vtable
size_t(* do_search)(struct interp_method *method, struct yac_interp_grid *grid, size_t *tgt_points, size_t count, struct yac_interp_weights *weights)
struct remote_point * data
struct weight_vector_3d * src_cell_gradient
struct supermesh_cell::@10 tgt
struct supermesh_cell::@10 src
struct weight_vector_data_3d * data
const_int_pointer num_vertices_per_cell
const_size_t_pointer cell_to_vertex_offsets
const_size_t_pointer cell_to_vertex
const_yac_int_pointer ids[3]
yac_const_coordinate_pointer vertex_coordinates
enum yac_edge_type * edge_type
double(* coordinates_xyz)[3]
void yac_quicksort_index_yac_int_size_t(yac_int *a, size_t n, size_t *idx)
void yac_quicksort_index_int_size_t(int *a, size_t n, size_t *idx)
static void yac_remove_duplicates_size_t(size_t *array, size_t *n)
void yac_quicksort_index_size_t_int(size_t *a, size_t n, int *idx)
#define YAC_ASSERT(exp, msg)
double(* yac_coordinate_pointer)[3]