(Eds)", Proceedings of the Computational Plasticity IV, CIMNE, Barcelona, 547-561. (1995), "Recent advances in the analysis and numerical simulation of strain localization in inelastic solids.
(2002), "On the use of embedded discontinuity elements with crack path continuity for mode-I and mixed-mode fracture", Eng. I think this should be clarified (perhaps on a case by case basis, gmsh is fine, invert order for Abaqus, etc) if we are to use imported meshes with JuAFEM. The Jacobi det is positive only in the clockwise case. In the meantime, I suspected that the problem is node ordering, so I made this example to reproduce it: using JuAFEM, TensorsĬ1 = Ĭ3 = Ĭoords = I will look into gmsh, and yes, that gear and it’s top and bottom surfaces should look like this. Thank you for your insights, they came just at the right moment, as right now I am looking into a linear elasticity problem on a 2d quadrilateral grid that came from Abaqus, and I ran into the same negative Jacobi determinants issue. It’s mathematically coherent and correct to do so, although you may need orientation elsewhere. That’s why I would recommend switching to absolute value of the volume for stiffness assembly in general. They should all be having the same orientation anyways, so it makes sense to be working with absolute value for stiffness assembly. By taking the absolute value of the n-submanifold you get the scalar volume of the element.
The volume of a full dimensional simplex is an n-submanifold, which for an embedded triangle in a higher dimension comes out as a quaternion for example. Since i need the volume to be a scalar and not a quaternion or other submanifold, it actually turns out to be correct to take the absolute value of the calculated volume. In my finite element assembly Adapode.jl, I actually fully switched to using the absolute value of volumes (for things like stiffness assembly) because for an embedded triangle inside a tetrahedral mesh, the volume comes out as a linear combinations of planes (basically isomorphic to a quaternion value). Unfortunately, I handled my unstructured grids by gmsh, therefore I’m not such a great help, but what I do know is, that you find here some references about the global-local indexing of 2d tets You can analyse this in detail with Debugger.jl Therefore, some elements are oscillating around 0 in terms of their det(J). Closely after, so one of the first elements, already results in a negative det(J) -1.7208e-6. Just try to understand the reference in the issue.Īnyhow, I tested the heat equation example and you have already in the first (and I think also in the second) element a det(J) of 1e-6, which is fairly close to zero.
I copy pasted this from 2 sources, so be aware of bugs in the code snippet above. In a test folder in your repo and people can then reproduce this easily. So something like grid = import_grid("gear.xml",output=:juafem_grid)Īddfaceset!(grid, "bottom", x -> norm(x) ≈ 0.025) Instead of testing it this way, I’d recommend to setup a test script, where you load your grid (“gears.xml”) and then try to export the dirichlet BC you want to set, see My guess would be (assuming I did everything right) is that the Fenics representation doesn’t observe the same face orientation rules as JuAFEM does, because every row in the connectivity matrix is in ascending order (while grid_generator creates ‘unordered’ tetrahedrons, probably with common faces having the same orientation in both tets). I don’t expect anyone to try and reproduce this, as that would be quite time consuming, but maybe someone has an idea just by reading this. If I try it with JuAFEM’s own grid = generate_grid(Tetrahedron, (10,10,10)) The issue is that trying to assemble the stiffness matrix triggers a det(J) is not positive positive error. This grid should be usable for the heat equation example almost directly, only by setting dim to 3, and handling the boundary conditions differently: ch = ConstraintHandler(dh) ĭbc_bottom = Dirichlet(:u, getfaceset(grid,"bottom"), (x,t) -> 0)ĭbc_top = Dirichlet(:u, getfaceset(grid,"top"), (x,t) -> 0) Adding them to the mesh above can be done by addfaceset!(grid, "bottom", x -> norm(x) ≈ 0.025)Īddfaceset!(grid, "top", x -> norm(x) ≈ 0.413226) So I can say grid = import_grid(path,output=:juafem_grid)Īnd it returns a Grid with the boundary matrix defined, but no facesets. only tetrahedral meshes), but should be sufficient for what I want to do, which is solving the heat equation on a gear mesh, as a proof of concept / learning example. I wanted to be able to use Fenics meshes, so I created a package to import them. Also apologies for not being able to create a copy-pastable MWE. Sorry in advance for posting this here, as this is most likely not a JuAFEM issue, but a problem with an external mesh.