For more complex boundaries, when creating submeshes, the lines at the top become quite complicated, but there are no problems with straight lines. How should I handle this?
from subprocess import *
import matplotlib.pyplot as plt
from dolfin import *
import os, shutil
import numpy as np
import time
import matplotlib
start = time.time()
parameters["form_compiler"]["quadrature_degree"] = 1
parameters["form_compiler"]["cpp_optimize"] = True
#网格
mesh = Mesh("brick.xml")
block = MeshFunction("size_t",mesh,"brick_physical_region.xml")
interface = MeshFunction("size_t",mesh,"brick_facet_region.xml")
submesh1 = MeshView.create(interface,8)
plot(submesh1)
plt.show()
In fact, I want to add an interface element between two bricks. In gmsh, I have defined a physical group for the interface, which is numbered 8, as shown in the sub-mesh of the figure below. However, when I output the undefined physical group numbers, the situation shown in the figure below occurred. The geometry I defined in gmsh is as follows:
As shown in the above diagram, the interface between the two bricks is defined as 8, and the rest of the areas are defaulted to 1. When importing the XML file and outputting the submesh, the physical groups displayed are 8 and 0, as shown above. There are no issues when outputting 8, but problems arise when outputting 0, resulting in a very messy appearance
Please add the actual Gmsh scripts to create your geometry to the post (so that people can reproduce your results).
What do you mean by this?
I would also like to see what happens if you output your submesh to xdmf, rather than plotting with matplotlib (as matplotlib plotting might produce weird results for «manifold» geometries.
When importing into the xdmf file, the results appear to be good. However, I am currently encountering another issue where I need to calculate the relative displacement at the interface to apply a cohesive force model. Could you please advise on how to apply this at the interface? Thank you very much for your answer
Without any code or a more detailed description, I do not know how to advise you on this. For instance, what do you imply by the “relative displacement”, is that:
An analytical expression
A finite element function resulting from solving a PDE
Thank you very much. I’ve learned a lot from your code. However, there’s one aspect I’m not very familiar with at the moment. Regarding the implementation of cohesion, what are the differences between globally applying DG (Discontinuous Galerkin) elements and using CG (Cohesive Zone) elements with the sub - grid method? In principle, what points should be noted when using DG elements? Or could you recommend some relevant articles for me to study?
First thing to note is that the DG formulation will always be more expensive as you have more degrees of freedom. Then it really depends on your application, do you know for sure that the crack will only occur at some interface because the surrounding materials are much more resistant or do not crack because they are ductile, then you can use CZM only on the interface. If you want to simulate possible crack propagation everywhere, I would recommend the DG approach. Note however that you will be very sensitive to mesh orientation. This is why phase field approaches tend to predominate nowadays.
Thank you for your reply. I currently want to simulate the failure behavior of composites. For the matrix, I plan to use the phase - field method to simulate crack propagation, while at the interface, I intend to insert interface elements to simulate interface failure, specifically as described in this paper. However, since I need to set up a function space for each brick to represent the relative displacement between every two bricks, it becomes very complicated. Therefore, I’m considering globally applying DG elements. What suggestions do you have?
The number of subdomains you have to consider is related to the maximum number of subdomains in “contact” with a given material. With the 4 color theorem, you can always find a pattern with 4 different subdomains in 2D, e;g. for the running bond pattern:
1 - 2 - 1 - 2 -
- 3 - 4 - 3 -
2 - 1 - 2 - 1 -
However, I understand that it may become a bit more complicated. DG everywhere is therefore simpler.
Thank you. However, I don’t quite understand the quadrature element at present. As shown in the following code, is it because the damage variable is in the UFL format and we need to obtain its vector, so we have to do it this way? But what does “quadrature” mean here?
Yes, since we decouple the displacement and damage problem, damage is represented as a Function living on the facet space. We would like to interpolate the UFL expression d_expr into this function space. Unfortunately, there is no interpolation mechanism yet for facet-based expressions, so we have to do it manually. This is done by assembling the vector corresponding to facet_interp which yields the wanted values by choosing a specific numerical integration (or quadrature) scheme as discussed in the tutorial.
Thank you for your advice. I now have a general understanding of how to solve this type of problem. I am using the Cohesive Zone Model only on the interface to handle complex interfaces, which might be complicated. If successful, I will open-source the code.
However, I currently have two questions. First, since I need to consider both tensile and compressive cases, there is an issue with the sign of the normal direction. In your code, you mentioned the problem regarding the value of effective_opening, but it hasn’t been considered here. If I need to take this into account, how should I approach it? Additionally, how is the positive direction of the normal vector specified? I am not familiar with this—is it based on the positive and negative signs we define?
Second, for the integral measure dINT, can subdomain_data choose MeshTags corresponding to facets? I found that this seems possible because defining subdomain_data in that way feels uncommon.
