underworldcode / UWGeodynamics

Underworld Geodynamics
Other
81 stars 32 forks source link

bug fix in _freesurface.py and add examples #250

Closed NengLu closed 3 years ago

NengLu commented 3 years ago
  1. fix a bug in _freesurface.py. TField should be set back to the initial mesh before advect the top nodes and remesh.
  2. add 3 examples to docs/examples: (1) Kaus2010_Rayleigh-Taylor_Instability (2) Crameri2012Case1_Relaxation (3) Crameri2012Case2_Rising_Plume

Reference

  1. Kaus, B. J., Mühlhaus, H., & May, D. A. (2010). A stabilization algorithm for geodynamic numerical simulations with a free surface. Physics of the Earth and Planetary Interiors, 181(1-2), 12-20.
  2. Crameri, F., Schmeling, H., Golabek, G. J., Duretz, T., Orendt, R., Buiter, S. J. H., ... & Tackley, P. J. (2012). A comparison of numerical surface topography calculations in geodynamic modelling: an evaluation of the ‘sticky air’method. Geophysical Journal International, 189(1), 38-54.
Results of Kaus2010_Rayleigh-Taylor_Instability at Time = 5 Ma (dt = 2.5 kyr) before after
image2000 image2000
rbeucher commented 3 years ago

So copying the mesh did not work?

NengLu commented 3 years ago

So copying the mesh did not work? mesh.data.copy() didn't work. and there are no function of mesh.copy() or copy(mesh),

rbeucher commented 3 years ago
from copy import copy

copy(mesh)
NengLu commented 3 years ago
from copy import copy

copy(mesh)

Humm, weird, that still didn‘t work. image Here are the mid of the RTI model mesh nodes coords, which should be around 0.500 while it deforms too much. image

rbeucher commented 3 years ago

Try with

from copy import deepcopy

NengLu commented 3 years ago

Try with

from copy import deepcopy

deepcopy can't work on model.mesh object, shows error: TypeError: can't pickle SwigPyObject objects

NengLu commented 3 years ago

Can we create the init_mesh in the _model.py below the def freeSurface, make it only works when the users using the freeSurf=True. It seems an easy way to build when the copy(mesh) didn't work, and more clear than the changes in commit https://github.com/underworldcode/UWGeodynamics/commit/43374a2c0bee06987784f34a0af29b157d763672

    def freeSurface(self, value):
        if value:
            self._init_mesh = 
            self._freeSurface = FreeSurfaceProcessor(self)
rbeucher commented 3 years ago

No because it goes against encapsulation. I would rather have the code as it is now.

NengLu commented 3 years ago

No because it goes against encapsulation. I would rather have the code as it is now.

Okay, thanks.