Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

rotation angles on corner grid points are wrong in cpld_gridgen #921

Closed
DeniseWorthen opened this issue Mar 15, 2024 · 0 comments · Fixed by #922
Closed

rotation angles on corner grid points are wrong in cpld_gridgen #921

DeniseWorthen opened this issue Mar 15, 2024 · 0 comments · Fixed by #922
Assignees
Labels
bug Something isn't working

Comments

@DeniseWorthen
Copy link
Contributor

DeniseWorthen commented Mar 15, 2024

Since MOM6 creates the model grid at runtime (including adjusting the land mask, if required), files required or used for the global applications for CICE and ATM must be created in a pre-processing step using only the MOM6 supergrid, topography and land mask files as input. This allows the mapped ocean mask (used for creating of the ATM ICs) and the CICE6 grid and mask files to be consistent with the run-time configuration of MOM6. A version of the cpld_gridgen utility has used to generate the required files even prior to being incorporated into ufs_utils, dating back at least 4 years when the authoritative repo was ufs-s2s.

For CICE, one of the required fields is the rotation angle for vector fields oriented I-J to be mapped to E-W and back. For CICE, this rotation angle is defined on the corner grid points, referred to as the Bu grid by MOM6. MOM6 itself calculates the rotation angle only on the center (Ct) grid points. In cpld_gridgen, in order to calculate the Bu angle, the Ct angle procedure from MOM6 was modified. This modification has been used since the earliest versions of cpld_gridgen.

This does not produce the correct rotation angles for CICE. This figure shows the Ct and (-1)*Bu angle in the current code

baseline

A slice through the middle of the Arctic shows the issue:
Screenshot 2024-03-15 at 1 58 43 PM

The Bu angle (black) and Ct angles should be of the same magnitude, since they are simply at slightly different geographic locations. But Bu is too small away from the tripole seam. For the 1/4deg tripole, the maximum angular error is ~20deg.

The fix is to calculate the Bu angle using an analogous procedure to that used internally by CICE to go from the Bu to the Ct angles. This procedure is shown as

  !-----------------------------------------------------------------
  ! Compute ANGLE on T-grid
  !-----------------------------------------------------------------
  if (trim(grid_type) == 'cpom_grid') then
     ANGLET(:,:,:) = ANGLE(:,:,:)
  else if (.not. (l_readCenter)) then
     ANGLET = c0


     !$OMP PARALLEL DO PRIVATE(iblk,i,j,ilo,ihi,jlo,jhi,this_block, &
     !$OMP                     angle_0,angle_w,angle_s,angle_sw)
     do iblk = 1, nblocks
        this_block = get_block(blocks_ice(iblk),iblk)
        ilo = this_block%ilo
        ihi = this_block%ihi
        jlo = this_block%jlo
        jhi = this_block%jhi


        do j = jlo, jhi
        do i = ilo, ihi
           angle_0  = ANGLE(i  ,j  ,iblk) !   w----0
           angle_w  = ANGLE(i-1,j  ,iblk) !   |    |
           angle_s  = ANGLE(i,  j-1,iblk) !   |    |
           angle_sw = ANGLE(i-1,j-1,iblk) !   sw---s
           ANGLET(i,j,iblk) = atan2(p25*(sin(angle_0)+ &
                                         sin(angle_w)+ &
                                         sin(angle_s)+ &
                                         sin(angle_sw)),&
                                    p25*(cos(angle_0)+ &
                                         cos(angle_w)+ &
                                         cos(angle_s)+ &
                                         cos(angle_sw)))
        enddo
        enddo
     enddo

This same procedure can be used to go from Ct to Bu with care for the indexing change and with knowing the Ct angle on the other side of the seam. Implementing this to calculate the Bu angles from the known MOM6 Ct angles produces the following

fix

Screenshot 2024-03-15 at 2 09 19 PM

In addition, once the Bu angles are known, the Ct angle that will be calculated by CICE (at run time) can again be calculated and that angle compared with the Ct angles produced by MOM6. This angchk has been implemented in a feature branch. The feature branch has been tested and shows that

  • the Ct angle produced by cpld_gridgen are O~10-6 or less different than the actual angle values produced by MOM6 at run time and written to the ocean history file. This has always been the case for the cpld_gridgen utility.
  • the Ct angle calculated using the Bu angles in cpld_gridgen (the angchk array) are O~10-8 different than the actual anglet values produced by CICE at runtime and written to the ice history file. This shows that the Bu angles calculated by cpld_gridgen are now consistent with the Ct angles from MOM6.
@GeorgeGayno-NOAA GeorgeGayno-NOAA added the bug Something isn't working label Mar 15, 2024
GeorgeGayno-NOAA pushed a commit that referenced this issue Apr 10, 2024
Remove find_anq routine and instead calculate Bu angles using 
the Ct angle on the opposite side of the seam.  Add angchk to 
verify CICE internally calculated Ct angle against both MOM6 and 
CICE model output.  Remove resetting of longitudes to 0:360.  Add  
tripole:tripole bilinear mapping for creation of CICE IC downscaling 
and add required Bu<->Ct mapping on tripole grid. Clean up log
messages, remove duplicate weight file and update documentation. 

Fixes #921.
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
bug Something isn't working
Projects
Development

Successfully merging a pull request may close this issue.

2 participants