Second-order Quadrilaterals (2D) meshes in standard Abaqus `.inp` format #2217

DanielDoehring · 2024-12-20T13:09:13Z

This implements second-order 2D quadrilateral p4est meshes constructed from a standard-format Abaqus .inp files.
See #1847 for this feature request.

Essentially, this PR is a bunch of tedious mesh file parsing to extract the relevant information. As in the existing implementation, we focus on meshes constructed e.g. from gmsh which are exported to .inp format.

For validation, I ran the simulation of a laminar, low Machnumber (essentially incompressible) simulation over the SD7003 airfoil.

The averaged drag & lift coefficient are

Average drag coefficient: 0.049819536195680685
Average lift coefficient: 0.3817639249290595

which are in excellent agreement with reference data (See this preprint)

github-actions · 2024-12-20T13:09:27Z

Review checklist

This checklist is meant to assist creators of PRs (to let them know what reviewers will typically look for) and reviewers (to guide them in a structured review process). Items do not need to be checked explicitly for a PR to be eligible for merging.

Purpose and scope

The PR has a single goal that is clear from the PR title and/or description.
All code changes represent a single set of modifications that logically belong together.
No more than 500 lines of code are changed or there is no obvious way to split the PR into multiple PRs.

Code quality

The code can be understood easily.
Newly introduced names for variables etc. are self-descriptive and consistent with existing naming conventions.
There are no redundancies that can be removed by simple modularization/refactoring.
There are no leftover debug statements or commented code sections.
The code adheres to our conventions and style guide, and to the Julia guidelines.

Documentation

New functions and types are documented with a docstring or top-level comment.
Relevant publications are referenced in docstrings (see example for formatting).
Inline comments are used to document longer or unusual code sections.
Comments describe intent ("why?") and not just functionality ("what?").
If the PR introduces a significant change or new feature, it is documented in NEWS.md with its PR number.

Testing

The PR passes all tests.
New or modified lines of code are covered by tests.
New or modified tests run in less then 10 seconds.

Performance

There are no type instabilities or memory allocations in performance-critical parts.
If the PR intent is to improve performance, before/after time measurements are posted in the PR.

Verification

The correctness of the code was verified using appropriate tests.
If new equations/methods are added, a convergence test has been run and the results
are posted in the PR.

Created with ❤️ by the Trixi.jl community.

codecov · 2024-12-20T14:04:11Z

Codecov Report

Attention: Patch coverage is 92.11823% with 16 lines in your changes missing coverage. Please review.

Project coverage is 90.05%. Comparing base (f09a707) to head (35100a4).

Files with missing lines	Patch %	Lines
src/meshes/p4est_mesh.jl	91.49%	16 Missing ⚠️

❗ There is a different number of reports uploaded between BASE (f09a707) and HEAD (35100a4). Click for more details.

HEAD has 2 uploads less than BASE

Flag BASE (f09a707) HEAD (35100a4)

unittests 26 24

Additional details and impacted files

@@            Coverage Diff             @@
##             main    #2217      +/-   ##
==========================================
- Coverage   96.37%   90.05%   -6.32%     
==========================================
  Files         486      487       +1     
  Lines       39186    39356     +170     
==========================================
- Hits        37764    35442    -2322     
- Misses       1422     3914    +2492

Flag	Coverage Δ
unittests	`90.05% <92.12%> (-6.32%)`	⬇️

Flags with carried forward coverage won't be shown. Click here to find out more.

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

sloede · 2024-12-20T16:35:11Z

This implements second-order 2D quadrilateral p4est meshes constructed from a standard-format Abaqus .inp files.

This is so awesome! Thanks a lot for this very cool feature 💪♥️

andrewwinters5000

Thanks for adding this feature! Overall the new file reading and constructors look good. I just left some questions about naming and future functionality.

NEWS.md

andrewwinters5000 · 2024-12-22T10:13:31Z

docs/src/meshes/p4est_mesh.md

@@ -386,6 +386,7 @@ transfinite map of the straight sided hexahedral element to find

 Also for a mesh in standard Abaqus format there are no qualitative changes when going from 2D to 3D.
 The most notable difference is that boundaries are formed in 3D by faces defined by four nodes while in 2D boundaries are edges consisting of two elements.
+Note that standard Abaqus also defines quadratic element types. In Trixi.jl, these higher-order elements are currently only supported in 2D, i.e., 8-node quadrilaterals.


Are only 8-node quads available? Or does this PR implement the infrastructure such that higher order with more nodes per element are (in theory) supported?

So at the moment yes, see lines 1931 to 1933 in p4est_mesh.jl. This should be easily generalizable, but I am actually not sure if there are elements with order higher than quadratic defined, see e.g. http://130.149.89.49:2080/v2016/books/usb/default.htm?startat=pt06ch28s01ael02.html

examples/p4est_2d_dgsem/elixir_navierstokes_SD7003airfoil.jl

andrewwinters5000 · 2024-12-22T10:26:12Z

src/meshes/p4est_mesh.jl


- `p4est_mesh_from_hohqmesh_abaqus`: High-order, curved boundary information created by
+- `p4est_mesh_from_hohqmesh_abaqus`: High-order, polygonial boundary information created by


I would prefer to say curved here. My main reason is that this is the word used in the HOHQMesh and other software (like NekMesh) documentation to discuss the high-order boundary information. I have never encountered anyone referring to it as polygonal. Also, "polygonal" here makes it sound like you are representing the boundary information with polygons, which is not what is happening.

src/meshes/p4est_mesh.jl

andrewwinters5000 · 2024-12-22T10:30:52Z

src/meshes/p4est_mesh.jl

+                            if occursin(quadratic_quads, current_element_type)
+                                # Print the first (element), second to fifth (vertices 1-4) indices to file.
+                                # For node order of quadratic (second-order) quads, see
+                                # http://130.149.89.49:2080/v2016/books/usb/default.htm?startat=pt06ch28s01ael02.html


Is this the ABAQUS user guide?

Yes, user guide/documentation. Unfortunately, I could not find the docs of more recent versions open accessible.

andrewwinters5000 · 2024-12-22T10:40:55Z

src/meshes/p4est_mesh.jl

+    if mesh_polydeg == 2
+        mesh_nnodes = mesh_polydeg + 1 # = 3
+        # Note: We ASSUME that the additional node between the end-vertices lies 
+        # on the center on that line, such that we can use Chebyshev-Gauss-Lobatto nodes!
+        # For polydeg = 2, we have the 3 nodes [-1, 0, 1] (within the reference element).


Interesting...so in the future for element boundaries higher than quadratic would this become mesh_polydeg >= 2?

This comment is also hyper-specific for the polydeg=2 case. In general the CGL nodes are nonuniform in [-1,1], include the endpoints, and will only include the point x=0 for even polydeg (i.e. an odd number of points).

Yeah so this comment is heavily tailored towards the quadratic case, which is as far as I now the highest degree supported by standard Abaqus (see also my answer above).

andrewwinters5000 · 2024-12-22T10:42:36Z

src/meshes/p4est_mesh.jl

+    if mesh_polydeg == 1
+        calc_tree_node_coordinates!(tree_node_coordinates, nodes, mapping,
+                                    vertices, tree_to_vertex)
+    else # mesh_polydeg = 2 => Nonlinearity in supplied mesh


okay, so this logic could readily handle higher order element orders (at least in 2D), right?

Yes, so more nodes here should be readily possible. I am not sure, however, if such elements are defined by the Abaqus standard.

andrewwinters5000 · 2024-12-22T10:49:21Z

src/meshes/p4est_mesh.jl

    n_trees = last(size(node_coordinates))
    nnodes = length(nodes)

    # Setup the starting file index to read in element indices and the additional
-    # curved boundary information provided by HOHQMesh.
+    # higher-order, polygonial boundary information provided by HOHQMesh.


Again, I would prefer to call the boundaries curved rather than polygonal. I already mentioned two reasons above. Another is that it over complicates the nomenclature. I think you want to use "polygonal" to clarify that the boundary is a closed, simply connected region with some number of sides (that are possibly individually curved). That is, some segments of the polygon can be straight-sided, attached to a spline and then attached to circular arc. This is why I think the catch all term of "curved" is better because it encompassed all these cases (as well as dimensions).

After all, one could say that a simple straight sided box domain has "polygonal boundary information" for the same reasons, but people would think I am weird to refer to a Carteisan box in this way :)

andrewwinters5000 · 2024-12-22T10:53:33Z

src/meshes/p4est_mesh.jl

+
+                    # The nodes for an Abaqus element are labeled following a closed path 
+                    # around the element:


Hmm, I imagine that figuring out this variable flipping will be significantly more cumbersome in 3D and for higher node orders.

Yeah, that is also why the 3D version is not shipped in this PR. I think one needs in that case to flip the faces also in some way which I have not yet figured out.

NEWS.md

src/meshes/p4est_mesh.jl

Co-authored-by: Andrew Winters <[email protected]>

DanielDoehring added 26 commits December 12, 2024 15:52

First steps towards quadratic/curved elements form standard Abaqus

e0e9aa4

idea

101448d

Notes

0b7d646

preproc abaqus

1e9af30

shorten

49d79c4

work on quad mesches

2c597a4

Seemingly working version 2nd order Abaqus 2D

70d46ca

Comments

7838c24

experiment with viz

b85e41b

test hohqmesh

82c506a

revert

c1f0d6c

comments

3bf1363

SD7003

fde8369

Merge branch 'main' into QuadraticElementsAbaqus

0248611

start 3d

5592a39

continue

1c51559

Continue

1ddc717

Reduce to 2D

0bbafa0

example

32a78f1

todo

3acd22a

Mention that boundary is in general only higher-order, not curved

0d0c07d

rename

393424a

test + example

41e5c31

docs

fc2d1f1

Comment

dcab80c

reset project toml

09e92c6

DanielDoehring requested a review from andrewwinters5000 December 20, 2024 13:09

DanielDoehring added 2 commits December 20, 2024 14:09

news

8b4173f

typo

b2bd8cb

Merge branch 'main' into QuadraticQuadsAbaqus

a07d25c

andrewwinters5000 requested changes Dec 22, 2024

View reviewed changes

DanielDoehring commented Dec 23, 2024

View reviewed changes

NEWS.md Outdated Show resolved Hide resolved

DanielDoehring commented Dec 23, 2024

View reviewed changes

src/meshes/p4est_mesh.jl Outdated Show resolved Hide resolved

DanielDoehring commented Dec 23, 2024

View reviewed changes

src/meshes/p4est_mesh.jl Outdated Show resolved Hide resolved

DanielDoehring and others added 5 commits December 23, 2024 09:01

Update src/meshes/p4est_mesh.jl

b00ea54

Update NEWS.md

523de2b

Update src/meshes/p4est_mesh.jl

1588a43

Co-authored-by: Andrew Winters <[email protected]>

Update examples/p4est_2d_dgsem/elixir_navierstokes_SD7003airfoil.jl

13f2994

Co-authored-by: Andrew Winters <[email protected]>

Update src/meshes/p4est_mesh.jl

35100a4

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Second-order Quadrilaterals (2D) meshes in standard Abaqus `.inp` format #2217

Second-order Quadrilaterals (2D) meshes in standard Abaqus `.inp` format #2217

DanielDoehring commented Dec 20, 2024

github-actions bot commented Dec 20, 2024

codecov bot commented Dec 20, 2024 •

edited

Loading

sloede commented Dec 20, 2024

andrewwinters5000 left a comment

andrewwinters5000 Dec 22, 2024

DanielDoehring Dec 24, 2024

andrewwinters5000 Dec 22, 2024

andrewwinters5000 Dec 22, 2024

DanielDoehring Dec 23, 2024

andrewwinters5000 Dec 22, 2024

DanielDoehring Dec 24, 2024

andrewwinters5000 Dec 22, 2024

DanielDoehring Dec 23, 2024

andrewwinters5000 Dec 22, 2024

andrewwinters5000 Dec 22, 2024

DanielDoehring Dec 23, 2024


		- `p4est_mesh_from_hohqmesh_abaqus`: High-order, curved boundary information created by
		- `p4est_mesh_from_hohqmesh_abaqus`: High-order, polygonial boundary information created by


		# The nodes for an Abaqus element are labeled following a closed path
		# around the element:

Second-order Quadrilaterals (2D) meshes in standard Abaqus .inp format #2217

Are you sure you want to change the base?

Second-order Quadrilaterals (2D) meshes in standard Abaqus .inp format #2217

Conversation

DanielDoehring commented Dec 20, 2024

github-actions bot commented Dec 20, 2024

Review checklist

Purpose and scope

Code quality

Documentation

Testing

Performance

Verification

codecov bot commented Dec 20, 2024 • edited Loading

Codecov Report

sloede commented Dec 20, 2024

andrewwinters5000 left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Second-order Quadrilaterals (2D) meshes in standard Abaqus `.inp` format #2217

Second-order Quadrilaterals (2D) meshes in standard Abaqus `.inp` format #2217

codecov bot commented Dec 20, 2024 •

edited

Loading