From 93381b65b01667ba120f6dcd0af4ec1771541624 Mon Sep 17 00:00:00 2001 From: Marshall Ward Date: Tue, 21 Jan 2025 17:39:04 -0500 Subject: [PATCH] Fix rotation in set_coupler_type_data `rotate_array` in `set_coupler_type_data` was trying to rotate an array to another of different size when idim and jdim are present. Some compilers seemed to let this through, but it raised a double-deallocation error in GCC. I'm guessing it's because the rotation was allocating to a new implicit array which was automatically deallocated. But I did not confirm this. This was modified to rotate onto a new array of the same size. The idim and jdim are passed to CT_set_data, which (hopefully) sorts out the indexing. --- src/framework/MOM_coupler_types.F90 | 19 ++++++++++++------- 1 file changed, 12 insertions(+), 7 deletions(-) diff --git a/src/framework/MOM_coupler_types.F90 b/src/framework/MOM_coupler_types.F90 index b931a2ddd2..b76912beb6 100644 --- a/src/framework/MOM_coupler_types.F90 +++ b/src/framework/MOM_coupler_types.F90 @@ -453,14 +453,12 @@ subroutine set_coupler_type_data(array_in, bc_index, var, solubility, scale_fact ! as array_in [A] integer :: subfield ! An integer indicating which field to set. integer :: q_turns ! The number of quarter turns through which array_in is rotated - integer :: is, ie, js, je, halo q_turns = 0 ; if (present(turns)) q_turns = modulo(turns, 4) subfield = ind_csurf if (present(solubility)) then ; if (solubility) subfield = ind_alpha ; endif if (present(field_index)) subfield = field_index - halo = 0 ; if (present(halo_size)) halo = halo_size ! The case with non-trivial grid rotation is complicated by the fact that the data fields ! in the coupler_2d_bc_type are never rotated, so they need to be handled separately. @@ -468,12 +466,19 @@ subroutine set_coupler_type_data(array_in, bc_index, var, solubility, scale_fact call CT_set_data(array_in, bc_index, subfield, var, & scale_factor=scale_factor, halo_size=halo_size, idim=idim, jdim=jdim) elseif (present(idim) .and. present(jdim)) then - ! Work only on the computational domain plus symmetric halos. - is = idim(2)-halo ; ie = idim(3)+halo ; js = jdim(2)-halo ; je = jdim(3)+halo - call allocate_rotated_array(array_in(is:ie,js:je), [1,1], -q_turns, array_unrot) + call allocate_rotated_array(array_in, [1,1], -q_turns, array_unrot) call rotate_array(array_in, -q_turns, array_unrot) - call CT_set_data(array_unrot, bc_index, subfield, var, & - scale_factor=scale_factor, halo_size=halo_size) + + if (modulo(q_turns, 2) /= 0) then + call CT_set_data(array_unrot, bc_index, subfield, var, & + idim=jdim, jdim=idim, & + scale_factor=scale_factor, halo_size=halo_size) + else + call CT_set_data(array_unrot, bc_index, subfield, var, & + idim=idim, jdim=jdim, & + scale_factor=scale_factor, halo_size=halo_size) + endif + deallocate(array_unrot) else call allocate_rotated_array(array_in, [1,1], -q_turns, array_unrot)