Quad Utilities: Enhanced Stability and New Documentation by mgharamti · Pull Request #1089 · NCAR/DART

mgharamti · 2026-03-28T02:15:30Z

Description:

This PR improves the numerical robustness of quadrilateral interpolation in quad_utils_mod and introduces new module documentation.

For the interpolation code, Cramer's rule (for solving the 3×3 system inquad_bilinear_interp) is now replaced with Gaussian Elimination with partial Pivoting (GEP). The new solver improves numerical stability, particularly for highly distorted or nearly degenerate quadrilateral cells. The computational cost remains constant (3×3 system) and negligible relative to the overall interpolation workflow.

I also added a new rst documentation page for quad_utils_mod. The documentation provides an overview of the interpolation procedure with detailed description of the supported grids, interpolation methods, and workflow.

Fixes issue

Fixes issue #1081

Types of changes

Bug fix (non-breaking change which fixes an issue)
New feature (non-breaking change which adds functionality)
Breaking change (fix or feature that would cause existing functionality to not work as expected)
Documentation update

Documentation changes needed?

My change requires a change to the documentation.
- I have updated the documentation accordingly.

Tests

I ran the developer tests for the quad_utils and also tested the new changes on small grids using different distorted quads. I got identical results between previous and refactored implementations for well-conditioned grids. For pathological quad configurations, the new solver demonstrates improved accuracy (orders of magnitude reduction in error).

Checklist for merging

Updated changelog entry
Documentation updated
Update conf.py

Checklist for release

Merge into main
Create release from the main branch with appropriate tag
Delete feature-branch

Testing Datasets

Dataset needed for testing available upon request
Dataset download instructions included
No dataset needed

Replaced Cramer's rule with Gaussian elimination with partial pivoting. This improves numerical stability especially for severely distored quadrilaterals.

hkershaw-brown · 2026-04-13T12:51:42Z

+                             'failed. It could be a singularity or a very tiny pivot'
+      write(string2, '(2(A, F10.4))') 'Obs location -> lon: ', lon, ', lat: ', lat  
+      write(string3, '(2(A, 4F10.4))') 'Quad corners -> x: ', x_corners, ', y: ', y_corners
+      call error_handler(E_MSG, 'quad_bilinear_interp', string1, &


This message will only print on task 0, so if other tasks have failures you will not see the message. To have the message print on any task:
E_ALLMSG

nice catch! I'll update this

hkershaw-brown · 2026-04-13T16:50:24Z

-! compute deter of m
-!d = deter3(m)
-
 ! Solve the matrix for b, c and d


@mgharamti Does this method need the do_rotate namelist option?

I ask because the code I think is doing a clockwise rotation (vs counterclockwise as the comment stays)
#833 (comment)
so if the rotation is not needed for stability with this new method, I think we should remove the rotation code.

Helen, in my tests rotations didn't add anything. I find the new solver to be stable without the need for rotations. Would that apply to the CICE examples is the issue you shared? probably.

The do_rotate true/false did not improve the out-of-bounds for the CICE

hmm, ok. I don't know if users actively use this rotation option. For legacy reasons, we could keep it. It's defaulted to .false. anyways.

ok I will put it on the todo list for the future to remove, and we leave it in for this pull request.

hkershaw-brown · 2026-04-13T17:59:20Z

@@ -0,0 +1,200 @@
+.. _quad_utils_mod:


Nice docs, I think it would be good to link to the quad_utils_mod from the "required model mod routines" doc

Definitely here:
https://docs.dart.ucar.edu/en/latest/guide/required-model-mod-routines.html#required-model-mod-routines

Maybe here?:
https://docs.dart.ucar.edu/en/latest/guide/advice-for-new-collaborators.html#interpolation

Something like:
"DART provides tools in the quad_utils_mod module for both grid search and interpolation on logically rectangular grids, that is, grids indexed like a regular (i,j) array but potentially curvilinear in physical space. These quad_utils_mod routines can be used within model_interpolate() to avoid having to implement interpolation routines from scratch"

Good idea, will do.

hkershaw-brown · 2026-04-13T18:00:47Z

+Mask Handling
+^^^^^^^^^^^^^
+
+An optional mask can be provided for fully irregular grids. If any of the four


This is just a curiosity question, why do only irregular grids except masks?

Helen, good question. To me, I think of the regular grid's interpolation as index-based. We don't need mask to define geometry; the grid itself is valid everywhere. So, masking here is handled outside (e.g., point_on_land, quad_on_land, ...).

On the other hand, interpolation in the irregular grids is geometry-based (physical space) in which we need the 4 corners to fit a surface. If one corner is land (or invalid), this will corrupt the interpolation. In essense, if any of the corners is masked we do not interpolate.

We could add masks for regular grid, but it adds overhead. It's essential in the irregular case for correctness.

Added references to the ``quad_utils_mod`` in 2 pages: guide/ required-model-mod-routines and advice-for-new-collaborators Also, changed message reporting when quad_interp fails to be done across all tasks.

hkershaw-brown · 2026-04-15T16:44:29Z

Question on the tests:

also tested the new changes on small grids using different distorted quads. I got identical results between previous and refactored implementations for well-conditioned grids.

Are these test refering to test_quad_reg_interp (well-conditioned grids) and test_quad_irreg_interp ( different distorted quads) or do you have other tests?

I ran the developer tests for the quad_utils

I'm not sure what the test is in the developer tests for quad_utils are testing, there's not a pass fail. So I just compared the .txt output from main and mgharamti:improve_quad_utils, (and plotted the results #1097)

The results are different on the last digit between main and this branch for test_quad_irreg_interp. I do not expect they would be bitwise between the 2 algorithms, because this is the point of the pull request.
Both this branch and main give 'large interp residual' messages.

[hkershaw:work] (improve_quad_utils) > ./test_quad_irreg_interp | grep 'large interp residual' | wc -l
17049

[hkershaw:work] (main) > ./test_quad_irreg_interp | grep 'large interp residual' | wc -l
15953

gfortran --version
GNU Fortran (MacPorts gcc13 13.3.0_2+stdlib_flag) 13.3.0
Copyright (C) 2023 Free Software Foundation, Inc.

mgharamti · 2026-04-15T20:30:36Z

Helen, I built my own distorted quad tests. I also ran the the ones in the developer tests but I didn't find them very informative.

hkershaw-brown · 2026-04-15T21:22:46Z

Helen, I built my own distorted quad tests. I also ran the the ones in the developer tests but I didn't find them very informative.

cool, I figured you had something cooler from the plots you showed at standup.

hkershaw-brown

Approved!

mgharamti added 2 commits March 27, 2026 19:51

Quad bilinear Solver

c5fdf5b

Replaced Cramer's rule with Gaussian elimination with partial pivoting. This improves numerical stability especially for severely distored quadrilaterals.

Added documentation for the quad utilities

78fa081

mgharamti added Documentation Improvements or additions to documentation quad_utils interpolation labels Mar 29, 2026

hkershaw-brown reviewed Apr 13, 2026

View reviewed changes

Improved docs

89658f3

Added references to the ``quad_utils_mod`` in 2 pages: guide/ required-model-mod-routines and advice-for-new-collaborators Also, changed message reporting when quad_interp fails to be done across all tasks.

hkershaw-brown mentioned this pull request Apr 15, 2026

Replace quad_bilinear_interp with inverse distance weighting in quad_utils. #887

Closed

15 tasks

Merge branch 'main' into improve_quad_utils

b717877

hkershaw-brown approved these changes Apr 22, 2026

View reviewed changes

hkershaw-brown added the release! bundle with next release label Apr 22, 2026

Conversation

mgharamti commented Mar 28, 2026

Description:

Fixes issue

Types of changes

Documentation changes needed?

Tests

Checklist for merging

Checklist for release

Testing Datasets

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

hkershaw-brown commented Apr 15, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

mgharamti commented Apr 15, 2026

Uh oh!

hkershaw-brown commented Apr 15, 2026

Uh oh!

hkershaw-brown left a comment

Choose a reason for hiding this comment

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

2 participants

hkershaw-brown commented Apr 15, 2026 •

edited

Loading