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

Problems with Mesh creation from info obtained from MeshGet() when MeshCreateDual() is involved #309

Closed
theurich opened this issue Oct 16, 2024 · 1 comment
Assignees
Labels
bug Something isn't working

Comments

@theurich
Copy link
Member

There are a few strange issues I am seeing when trying to create a Mesh from info I pull from an existing Mesh via ESMF_MeshGet().

Background, and what works as expected

I start with a very basic E4P4 cubed sphere mesh generated by NEPTUNE for debugging. E4P4 means that each of the six faces of the cube hold 4x4 spectral elements in the horizontal, each of which will have 4 nodes of the basis function. The ESMF Mesh is constructed on those nodal points, leading to total of 6 x (4 x 4) x (4 x 4) + 2 = 1538 nodes in the Mesh. I distribute across 2 PETs which leads to manageable situation to look at visually (there will always be two pictures show, PET0 and PET1). The original Mesh looks as expected:

image
image

I use the following sequence of calls to pull out info from the original mesh, and construct a duplicate meshNew:

    call ESMF_MeshGet(mesh, nodeCount=nCount, elementCount=eCount, &
      elementConnCount=eConnCount, coordSys=coordSysOrig, &
      nodalDistgrid=nDistgrid, elementDistgrid=eDistgrid, &
      parametricDim=pdim, spatialDim=sdim, rc=rc); CHECKRC
write(0,*) "pdim=", pdim, "sdim=", sdim
    allocate(nodeIds(nCount), nodeCoords(nCount*2), nodeOwners(nCount))
    allocate(elementIds(eCount), elementTypes(eCount), elementConn(eConnCount))
    call ESMF_MeshGet(mesh, nodeIds=nodeIds, nodeCoords=nodeCoords, &
      nodeOwners=nodeOwners, elementIds=elementIds, elementTypes=elementTypes, &
      elementConn=elementConn, rc=rc); CHECKRC
    meshNew = ESMF_MeshCreate(parametricDim=pdim, spatialDim=sdim, &
      coordSys=coordSys, rc=rc) ; CHECKRC
    call ESMF_MeshAddNodes(meshNew, nodeIds=nodeIds, nodeCoords=nodeCoords, &
      nodeOwners=nodeOwners, nodalDistGrid=nDistgrid, rc=rc) ; CHECKRC
    call ESMF_MeshAddElements(meshNew, &
      elementIds=elementIds, elementTypes=elementTypes, elementConn=elementConn, &
      elementDistgrid=eDistGrid, rc = rc) ; CHECKRC

Doing this on the original NEPTUNE E4P4 mesh results in the expected duplication:

image
image

So far so good....

First odd behavior

The first strange thing I encounter is when I create the dual mesh for the E4P4 mesh. In order to do so I add element coords when I create the E4P4 mesh that I calculate as "center of mass" from the surrounding nodes. The mesh produced by ESMF_MeshCreateDual() looks acceptable in terms of shape and topology. The nodes shown are essentially the center of mass coords I filled into the original E4P4 mesh. However, there seems to be a distribution issue that can be seen at the interface between the PET0 and PET1 plots. The interface elements are replicated between both PETs. That is the first issue I do not really understand!

image

image

Second surprise - a pancake!

The second surprise I ran into is when I use the same code shown above to pull info from the dual mesh and try to construct a duplication of it. It now leads to strange flat structure:

image
image

It almost looks to me as if the coords come out in radians and not degree, although I checked and the coordSys returned from MeshGet() is ESMF_COORDSYS_SPH_DEG. Also, why would this suddenly be different than when I used the duplication code on the original E4P4 mesh? I use the same sequence of call, so if coords came out correctly the first time, why is it not working on the dual mesh?

Explicit redist fixes first issue, but not second one

In order to address the first issue of distribution of the dual mesh, I send it through ESMF_MeshCreate(fromMesh):

      meshDualRedist = ESMF_MeshCreate(meshDual, &
        nodalDistgrid=eDistgrid, elementDistgrid=nDistgrid, rc=rc); CHECKRC

where I eDistgrid is the elementDistgrid and nDistgrid is the nodalDistgrid from the original E4P4 mesh.

With that the distribution issue seems to go away, i.e. I do not see element duplication between PETs:

image
image

However, again using the same MeshGet/MeshCreate sequence to create a duplicate of this redistributed dual mesh still shows the same 2nd surprise of the flat pancake mesh:

image
image

Of course the replicated elements at the interface between the PETs is also gone here, as would be expected, but still it is a flat structure.

Attempt to fix the pancake issue

Applying a forced rad -> deg conversion on the nodeCoords array pulled out with MeshGet() for the dual mesh duplication, I do get a bit of the expected curvature back into the duplicated mesh. However, the two PETs are not covering the full sphere. Something is definitely much more messed up than just a rad vs deg issue. And also still, the original E4P4 mesh duplication via MeshGet() information works just fine, so why is the dual mesh affected?

image
image

@theurich theurich added the bug Something isn't working label Oct 16, 2024
@theurich
Copy link
Member Author

The pancake issue has been resolved on branch https://github.com/esmf-org/esmf/tree/fix_meshdual. The other issue of distribution after Mesh_CreateDual() will be addressed under #308.

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
None yet
Development

No branches or pull requests

2 participants