From baded5efea72c9d61caa51f5d2221366e2f9410c Mon Sep 17 00:00:00 2001 From: Patrick Scholz Date: Fri, 26 Jun 2026 14:33:05 +0200 Subject: [PATCH] Fix mEVP/aEVP: zero stress tensor before thin-ice cycle to prevent stale stress application stress_tensor_m and stress_tensor_a skip elements with low ice mass (msum <= 0.01) but leave sigma11/22/12 unchanged. stress2rhs_m uses a different skip condition (concentration < 0.01) and may not skip the same element. Stale stresses from when the ice was thicker are then applied to the momentum equation, causing spurious acceleration of nearly-vanished ice during summer decay. Explicitly zero the stress components before cycling so stress2rhs_m always applies zero contribution for elements below the mass threshold. --- src/ice_maEVP.F90 | 14 ++++++++++++-- 1 file changed, 12 insertions(+), 2 deletions(-) diff --git a/src/ice_maEVP.F90 b/src/ice_maEVP.F90 index aec61801c..d08a7fb1e 100644 --- a/src/ice_maEVP.F90 +++ b/src/ice_maEVP.F90 @@ -138,7 +138,12 @@ subroutine stress_tensor_m(ice, partit, mesh) if (ulevels(elem) > 1) cycle msum=sum(m_ice(elnodes))*val3 - if(msum<=0.01_WP) cycle !DS + if(msum<=0.01_WP) then + sigma11(elem) = 0.0_WP + sigma22(elem) = 0.0_WP + sigma12(elem) = 0.0_WP + cycle + end if asum=sum(a_ice(elnodes))*val3 dx=gradient_sca(1:3,elem) @@ -1061,7 +1066,12 @@ subroutine stress_tensor_a(ice, partit, mesh) elnodes=elem2D_nodes(:,elem) msum=sum(m_ice(elnodes))*val3 - if(msum<=0.01_WP) cycle !DS + if(msum<=0.01_WP) then + sigma11(elem) = 0.0_WP + sigma22(elem) = 0.0_WP + sigma12(elem) = 0.0_WP + cycle + end if asum=sum(a_ice(elnodes))*val3 dx=gradient_sca(1:3,elem)