Skip to content

Instantly share code, notes, and snippets.

@Jason-S-Ross
Last active August 5, 2021 14:05
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save Jason-S-Ross/11336ddb94f59cf8377e1ffe09df1db5 to your computer and use it in GitHub Desktop.
Save Jason-S-Ross/11336ddb94f59cf8377e1ffe09df1db5 to your computer and use it in GitHub Desktop.
Question on Tensor Calculus for Sympy User Group

I'm trying to use Sympy to help solve problems in solid mechanics, where I have a lot of tensor differential equations, expressed in arbitrary curvilinear coordinates in euclidean space.

Here's an example problem:

The compatibility of displacements requires that Show that the following is displacement field satisfies the compatibility requirement and compute the components of stress in spherical coordinates: where is a scalar field that is harmonic, i.e. and the stress is given by the following equation:

This creates a few challenges:

  • How do we manipulate differential equations involving covariant derivatives like symbolically, to show that a given equation like satisfies the differential equation without abandoning the index notation?
  • How do we take a known solution like and compute the corresponding stress in an arbitrary coordinate system like spherical coordinates? This involves computing the metric tensor for the coordinate system, computing Christoffel symbols, and evaluating covariant derivatives, all of which are possible, but the Tensor api doesn't seem to support.
  • How do we operate on tensor expressions involving the metric tensor itself, when replace_with_arrays relies on the metric being a property of the TensorIndexType? For instance, consider the following:
    from sympy import *
    from sympy.tensor.tensor import TensorHead, TensorSymmetry, TensorIndexType, tensor_indices
    from sympy.abc import eta
    euclid = TensorIndexType('Euclid', dummy_name='R')
    i, j, r, s = tensor_indices('i j r s', euclid)
    v = TensorHead('v', [euclid])
    # Hack to handle the first covariant derivative of v; populate with values "by hand"
    v1 = TensorHead('v', [euclid, euclid])
    # Hack to handle the second covariant derivative of v; populate with values "by hand"
    v2 = TensorHead('v', [euclid, euclid, euclid])
    g = TensorHead('g', [euclid, euclid], TensorSymmetry.fully_symmetric(2))
    expr = g(j, s) * v1(i, -s) + g(i, r) * v1(j, -r) + 2 * eta / (1 - 2 * eta) * g(i, j) * v1(r, -r)

How do we evaluate this expression in spherical coordinates? g is the metric tensor, but the substitution we give in the replacement dictionary for euclid is what's actually used to raise and lower indices.

  • How do we compute Christoffel symbols for tensor expressions? They're part of the diffgeom module, but they're not represented there in a way that's compatible with the tensor module.

In light of these I'm forced to work with piecemeal solutions that require finicky expression manipulation using combinations of tensorproduct and tensorcontraction and lose a lot of the expressiveness of tensor notation in the process. Can anyone provide some pointers on working with tensor differential equations in sympy?

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment