Maybe we should stop using planar triangles for fault modeling.

Let’s consider a model problem where the goal is to calculate the stress on a fault due to some slip on that fault. We also include a topographic free surface. How should we model this surface?

Should we use planar triangles to model a curved fault? What happens at the junction between two triangles? If we do the natural fault thing and require the slip to be tangential to the suface, the slip vector has a step function at the junction between the two triangles. This must be true because the two triangles have different normal vectors and tangent vectors. On the other hand, if we force slip to be continuous across triangles, we have to allow non-tangential slip because at an edge where two triangles meet the slip vector has to be identical but the tangent vectors are different.

That means that we’re stuck between a rock and a hard place and have to choose:

So, planar triangles probably aren’t a sustainable or “correct” way to represent a fault surface unless that surface is perfectly planar or we’re zoomed in far enough that the underlying microscale behavior is relevant. Fundamentally the problem is that our fault surface representation is $C^0$, meaning that the function itself is continuous but the derivatives are not. This problem propagates through the whole model up in other ways too. Suppose you choose to have discontinuous, but tangential, slip. Then, there will be stress singularities at every element junction!

Should we use a smoother fault surface representation? Probably. A $C^1$ fault representation where the normal vector is continuous between adjacent elements solves the problem described in the previous paragraph.

Let’s back up and consider how to mesh our model earthquake problem from the perspective of truth, accuracy and robustness:

Truth: The underlying true surface of the Earth is fractal and has sharp corners at basically any zoom level. The fault surfaces are probably similarly fractal and non-smooth.

Accuracy: There is accurate satellite data on the location of the Earth’s surface down to about a a few meters resolution. And for faults, depending on where in the world, the uncertainty in fault location is on the order of 10-10000m. So, only a few digits of accuracy could possibly be relevant in these models since there’s already a high level of error introduced by the mere location uncertainty in the faults, much less all the other uncertainties that I haven’t discussed here. The main goal, at least at this point in the evolution of the earthquake science field, is in terms of qualitative behavior and questions on the scale of 1-50% or even about the sign of an effect (positive or negative).

Robustness: This is where smoothness gets really important. From my experience, having smooth boundaries makes so much about integral equation methods simpler. Especially for these faulting/crack problems where hypersingular terms are common.

A different approach

My tentative approach that I’ve been exploring is:

Thanks to Andreas Klöckner, Brendan Meade and Andrew Bradley for discussions that led to this post. If you want to chat about these kind of topics, please get in touch!

Posts BIEBook