[BUG] - mesh.X.coords does not display the updated coordinates after mesh deforming

[NengLu]

After deforming the mesh, mesh.X.coords does not display the updated coordinates.
The coordinates in mesh.X.coords do not match those obtained from mesh.dm.

Results from the latest development branch:

Results from some version of the old development branch:

The bug could be tested by:

import underworld3 as uw
import numpy as np
import sympy

xmin, xmax = -0.5, 0.5
ymin, ymax = -1,0
xres, yres = 50,50
mesh = uw.meshing.StructuredQuadBox(elementRes=(int(xres), int(yres)), minCoords=(xmin, ymin), maxCoords=(xmax, ymax))

deform_fn = 0.1 * sympy.cos(2.*sympy.pi*(mesh.X[0])/1.0)
Dz = uw.discretisation.MeshVariable("Dz", mesh, 1, degree=1)
diffuser = uw.systems.Poisson(mesh, Dz)
diffuser.constitutive_model = uw.constitutive_models.DiffusionModel
diffuser.constitutive_model.Parameters.diffusivity = 1.
diffuser.add_essential_bc((deform_fn,), "Top")
diffuser.add_essential_bc((0.,),"Bottom")
diffuser.solve()
displacement = np.zeros((mesh.X.coords.shape[0],mesh.dim))
displacement[:,-1] = uw.function.evaluate(Dz.sym[0], mesh.X.coords)[:,0,0]
mesh._deform_mesh(mesh.X.coords + displacement)

import matplotlib.pyplot as plt
from underworld3.cython.petsc_discretisation import petsc_dm_find_labeled_points_local
topwall = petsc_dm_find_labeled_points_local(mesh.dm,"Top") 
coords = mesh.X.coords[topwall]
index =  np.argsort(coords[:,0])
xx = coords[index,0]
yy = coords[index,-1]
plt.plot(xx,yy,label=r'coords of topwall from mesh.X.coords')
coord_vec = mesh.dm.getCoordinatesLocal()
coords = coord_vec.array.reshape(-1, mesh.cdim)[topwall]
index =  np.argsort(coords[:,0])
xx = coords[index,0]
yy = coords[index,-1]
plt.plot(xx,yy,label=r'coords of topwall from mesh.dm.getCoordinatesLocal')
plt.legend(loc = 'upper left',prop = {'size':8})
plt.savefig('tem.png',dpi=300,bbox_inches='tight')

[NengLu]

Thanks, Louis. The fix works. I am closing it.

[NengLu]

I need to reopen this issue because the latest other merge introduced a related problem (if only merging this one would not cause the new problem).
After another mesh deformation, the mesh variable wasn't updated by the Stokes solver.

Could be tested by:

import underworld3 as uw
import numpy as np
import sympy
import matplotlib.pyplot as plt
from underworld3.cython.petsc_discretisation import petsc_dm_find_labeled_points_local

xmin, xmax = -0.5, 0.5
ymin, ymax = -1,0
xres, yres = 50,50
visc_fn = 4e-5
dt_set = 5e-5

mesh = uw.meshing.StructuredQuadBox(elementRes=(int(xres), int(yres)), minCoords=(xmin, ymin), maxCoords=(xmax, ymax))
v = uw.discretisation.MeshVariable("V", mesh, mesh.dim, degree=1,continuous=True)
p = uw.discretisation.MeshVariable("P", mesh, 1, degree=0,continuous=False)
Mxy = uw.discretisation.MeshVariable("Mxy", mesh, vtype=uw.VarType.SCALAR, degree=1, continuous=True, varsymbol=r"{\cal{Mxy}}")
 
deform_fn = 0.02 * sympy.cos(2.*sympy.pi*(mesh.X[0])/1.0)
Dz = uw.discretisation.MeshVariable("Dz", mesh, 1, degree=1)
diffuser = uw.systems.Poisson(mesh, Dz)
diffuser.constitutive_model = uw.constitutive_models.DiffusionModel
diffuser.constitutive_model.Parameters.diffusivity = 1.
diffuser.add_essential_bc((deform_fn,), "Top")
diffuser.add_essential_bc((0.,),"Bottom")
diffuser.solve()
displacement = np.zeros((mesh.X.coords.shape[0],mesh.dim))
displacement[:,-1] = uw.function.evaluate(Dz.sym[0], mesh.X.coords)[:,0,0]
mesh._deform_mesh(mesh.X.coords + displacement)

stokes = uw.systems.Stokes(mesh, velocityField=v, pressureField=p)
stokes.constitutive_model = uw.constitutive_models.ViscousFlowModel
stokes.bodyforce = sympy.Matrix([0, -1.0])
stokes.constitutive_model.Parameters.shear_viscosity_0 = visc_fn
stokes.saddle_preconditioner = 1.0 / stokes.constitutive_model.Parameters.shear_viscosity_0
stokes.add_essential_bc((0.0,None), "Left")
stokes.add_essential_bc((0.0,None), "Right")
stokes.add_essential_bc((0.0,0.0), "Bottom")
stokes.tolerance = 1.0e-6
stokes.petsc_options["ksp_rtol"] = 1.0e-6
stokes.petsc_options["ksp_atol"] = 1.0e-6
stokes.petsc_options["snes_converged_reason"] = None
stokes.petsc_options["snes_monitor_short"] = None

stokes.solve(zero_init_guess=False)
topwall = petsc_dm_find_labeled_points_local(mesh.dm,"Top") 
coords = mesh.X.coords[topwall]
index =  np.argsort(coords[:,0])
xx = coords[index,0]
yy = coords[index,-1]
vel = uw.function.evaluate(v,coords)[:,0,:].copy()
print(vel.max())
print(coords[:,-1].max())
print(stokes.mesh.X.coords[topwall,-1].max())
print(v.coords[topwall,-1].max())

deform_fn = 0.05 * sympy.cos(2.*sympy.pi*(mesh.X[0])/1.0)
Dz2 = uw.discretisation.MeshVariable("Dz2", mesh, 1, degree=1)
diffuser2 = uw.systems.Poisson(mesh, Dz2)
diffuser2.constitutive_model = uw.constitutive_models.DiffusionModel
diffuser2.constitutive_model.Parameters.diffusivity = 1.
diffuser2.add_essential_bc((deform_fn,), "Top")
diffuser2.add_essential_bc((0.,),"Bottom")
diffuser2.solve()
displacement = np.zeros((mesh.X.coords.shape[0],mesh.dim))
displacement[:,-1] = uw.function.evaluate(Dz2.sym[0], mesh.X.coords)[:,0,0]
mesh._deform_mesh(mesh.X.coords + displacement)

stokes.solve(zero_init_guess=False)
coords = mesh.X.coords[topwall]
index =  np.argsort(coords[:,0])
xx = coords[index,0]
yy = coords[index,-1]
vel = uw.function.evaluate(v,coords)[:,0,:].copy()
print(vel.max())
print(coords[:,-1].max())
print(stokes.mesh.X.coords[topwall,-1].max())
print(v.coords[topwall,-1].max())

The output is:

Structured box element resolution 50 50
  0 SNES Function norm 0.0197597
  Nonlinear Solver_13_ solve did not converge due to DIVERGED_MAX_IT iterations 0
  0 SNES Function norm 0.0197597
  1 SNES Function norm 1.17214e-05
  2 SNES Function norm < 1.e-11
  Nonlinear Solver_13_ solve converged due to CONVERGED_FNORM_ABS iterations 2
39.50464698012091
0.02
0.02
0.02
  0 SNES Function norm 0.0198329
  Nonlinear Solver_13_ solve did not converge due to DIVERGED_MAX_IT iterations 0
  0 SNES Function norm 0.149641
  1 SNES Function norm 3.47777e-07
  2 SNES Function norm < 1.e-11
  Nonlinear Solver_13_ solve converged due to CONVERGED_FNORM_ABS iterations 2
39.50464698012091
0.07
0.07
0.07

The max velocity on the topwall should be 131.91088367986202, while it keeps the same value as 39.50464698012091

[lmoresi]

@NengLu — proposed fix in PR #140, ready for review.

Root cause (re-opened portion)

Two holes combined, only visible after PR #127 ("Trust JIT cache: skip DM rebuild on constant-only parameter changes", merged 2026-04-21):

1. `_deform_mesh` never marked attached solvers `is_setup = False`. Pre-Trust JIT cache: skip DM rebuild on constant-only parameter changes #127 this was harmless because the solver rebuilt its PETSc DM on every `solve()` anyway. Post-Trust JIT cache: skip DM rebuild on constant-only parameter changes #127 the cached DM persists across solves, so coord changes have to be signalled explicitly. `mesh.adapt()` and submesh re-extraction already did this; `_deform_mesh` didn't.

2. `_build()` Fast Path 1 (constants-only update) short-circuited on a matching JIT cache key without checking `_needs_dm_rebuild`. Fast Path 2 already checks it; Fast Path 1 was missing the same guard. So even after fix (1) sets `is_setup = False`, Fast Path 1 would fire and keep the stale DM.

Symptom: `F(v_prev)` on the cached pre-deform DM ≈ 0 → SNES reports "converged in 0 iterations" → velocity field never updates. That's why your second deform kept `vel.max = 39.50` instead of scaling to 131.91.

Verification on your reproducer

	First deform	Second deform
Pre-fix	`vel.max = 39.50`	`vel.max = 39.50` (stuck)
Post-fix	`vel.max = 39.50`	`vel.max = 131.91` (matches the analytical value you quoted)

Ratio 3.34×, matching the expected linear scaling.

Also added `tests/test_0820_deform_mesh_solver_rebuild_regression.py` — a level_1 regression distilled from your reproducer — so this doesn't silently regress again.

Next

Would you be able to pull PR #140 (branch `bugfix/deform-mesh-solver-invalidate-122`) and re-run your original script? If the post-second-deform velocity behaves as you expect, I'll merge and close the issue.

Underworld development team with AI support from Claude Code

[NengLu]

Thanks, Louis. The quick test (using branch 'bugfix/deform-mesh-solver-invalidate-122') showed that the fix resolves the issue, and the old ALE and ALEIB codes are working again. The model test results are comparable to the analytical solution. It's okay to close the issue and proceed with the merge.

I am now moving forward with the ALEIB implementation using the 'submesh'.

[lmoresi]

Fixed in #140 (merged to development at 09d8b5a). NengLu confirmed working on the bugfix branch. Closing manually — GitHub auto-close didn't fire because the PR merged to development not main.