Where communities thrive

  • Join over 1.5M+ people
  • Join over 100K+ communities
  • Free without limits
  • Create your own community
Repo info
  • 05:10
    mohajain commented #22353
  • 05:00
    smichr commented #22090
  • 04:39
    auntyellow opened #22361
  • 03:02
    github-actions[bot] commented #22211
  • 02:11
    eagleoflqj commented #22169
  • 02:10
    eagleoflqj commented #22169
  • 01:49
    sympy-bot commented #22211
  • 01:49
    eagleoflqj synchronize #22211
  • Oct 24 22:10
    ThePauliPrinciple commented #22360
  • Oct 24 22:05
    ThePauliPrinciple commented #22360
  • Oct 24 22:03
    ThePauliPrinciple commented #22360
  • Oct 24 21:55
    ThePauliPrinciple commented #22360
  • Oct 24 20:15
    jkr26 starred sympy/sympy
  • Oct 24 18:41
    ayushk7102 commented #22275
  • Oct 24 18:40
    smichr closed #22342
  • Oct 24 18:40
    smichr commented #22342
  • Oct 24 18:33
    ThePauliPrinciple commented #22351
  • Oct 24 16:32
    qiayuanliao starred sympy/sympy
  • Oct 24 15:25
    Albuquerque-David starred sympy/sympy
  • Oct 24 13:30
    yasserotiefy starred sympy/sympy
Jason Ross

I'm having trouble learning the diffgeom module. I'm trying to compute covariant derivatives in spherical coordinates and I feel like I'm missing something.

Here's what I'm working with:

from sympy import *
from sympy.abc import x, y, z, r, theta, phi
from sympy.diffgeom import CoordSystem, Manifold, Patch, TensorProduct as TP

R3 = Manifold("R3", 3)
S = Patch("S", R3)
relations = {
    ("Car3D", "Sph") : Lambda(
        (x, y, z),
            sqrt(x**2 + y**2 + z**2), 
            acos(z / sqrt(x**2 + y**2 + z**2)),
            atan2(y, x),
    ("Sph", "Car3D"): Lambda(
        (r, theta, phi),
        Matrix([r * sin(theta) * cos(phi), r * sin(theta) * sin(phi), r * cos(theta)])
Car3D = CoordSystem("Car3D", S, [x, y, z], relations)
Sph = CoordSystem("Sph", S, [r, theta, phi], relations)

I can compute the metric tensor in spherical coordinates by doing the following:

g = (Car3D.jacobian(Sph) @ Car3D.jacobian(Sph).T).simplify()

However, this computes the metric tensor with respect to $x, y, z$, while I'd like to write the metric tensor with respect to $r, \theta, \phi$:
[10001x2+y2+z20001x2+y2]\left[\begin{matrix}1 & 0 & 0\\0 & \frac{1}{x^{2} + y^{2} + z^{2}} & 0\\0 & 0 & \frac{1}{x^{2} + y^{2}}\end{matrix}\right]
I can do a substitution, but it's kind of awkward:

x, y, z = Car3D.symbols
r, theta, phi = Sph.symbols
forward_subs = {
    x: r * sin(theta) * cos(phi),
    y: r * sin(theta) * sin(phi),
    z: r * cos(theta)
g_sph = g.subs(forward_subs).simplify()

[10001r20001r2sin2(θ)]\left[\begin{matrix}1 & 0 & 0\\0 & \frac{1}{r^{2}} & 0\\0 & 0 & \frac{1}{r^{2} \sin^{2}{\left(\theta \right)}}\end{matrix}\right]
This gives me the expected metric tensor as a matrix but there's a lot of extra effort:

  • I have to define the forward and backward transformations BEFORE I declare forward_subs since the CoordinateSymbol "x" does not exist until I create Car3D, and I have to provide the forward and backward transformations to create Car3D.
  • Converting from one coordinate system to another uses a dictionary that duplicates all the information already defined in the coordinate systems
  • Computing a metric tensor by taking the Jacobian times its transpose seems like extra effort
  • I still haven't got a tensor expression for the metric tensor

I can convert the metric tensor Matrix to a tensor expression as follows:

g_sph_tens = 0
for i, e_i in enumerate(Car3D.base_oneforms()):
    for j, e_j in enumerate(Car3D.base_oneforms()):
        g_sph_tens += g_sph[i, j] * TP(e_i, e_j)

However, computing the Christoffel symbol for g_sph_tens using metric_to_Christoffel_2nd gives all zeros, which is incorrect.

Can somebody help me out with this? It seems like I'm missing some fundamentals.

1 reply
Matías Senger

Hi, I want to learn SymPy to do some calculations. I am practicing with the electromagnetic Lagrangian term $F_{\mu\nu}F^{\mu\nu}$. I have defined $F$ in the way that is explained in the documentation but cannot figure out how to go from there to $E^2 + B^2$. How am I supposed to do this?

Up to now I have this:

Ex, Ey, Ez, Bx, By, Bz = sp.symbols('E_x E_y E_z B_x B_y B_z')
c = sp.symbols('c', positive=True)

F = spT.TensorHead('F', [Lorentz, Lorentz], spT.TensorSymmetry.fully_symmetric(-2))

repl = {Lorentz: sp.diag(1, -1, -1, -1)}
repl.update({F(-mu,-nu): [
    [0, Ex/c, Ey/c, Ez/c],
    [-Ex/c, 0, -Bz, By],
    [-Ey/c, Bz, 0, -Bx],
    [-Ez/c, -By, Bx, 0]]}

and I have tried (F(mu,nu)*F(-mu,-nu)).replace_with_arrays(repl, [mu,nu]) and also F(mu,nu).replace_with_arrays(repl, [mu,nu])*F(-mu,-nu).replace_with_arrays(repl, [mu,nu]).

6 replies
Priyanshu Agarwal
Hey, I have some huge equations but they are getting printed incomplete and half when I use sympy preview for pdf (see the pdf here) , could you please help me with this
sympy.preview(expr, output="pdf", viewer="file", filename="DFN_equations.pdf", dvioptions=['-D', '1100'], euler=False)
Jason Ross

I wanted to expand on my question from earlier since it may be too specific.
I am working through some equations from the theory of elasticity.

The differential equation I am trying to work with is the equation
for compatibility of small displacements in a strained isotropic elastic body.
This is written as
(12η)vijj+vjji=0(1 - 2\eta) v^i|^j_j + v^j|^i_j = 0
represents the covariant derivative of f with respect to j, and
represents the covariant derivative of f with respect to j
with the j index raised (this convention is from Green and Zerna -
I don't know if this is a standard notation).

We can use solutions of the above differential equation to compute
the stress tensor as follows:

τijμ=gjsvis+girvjr+2η12ηgijvrr\frac{\tau^{ij}}{\mu} = g^{js} v^i|_s + g^{ir} v^j|_r + \frac{2 \eta}{1 - 2 \eta} g^{ij} v^r|_r

One of the solutions to the displacement equation is as follows:

vi=Fiv_i = F|_i

where F is a harmonic function.

What's a good way to

  1. Verify that this satisfies the displacement equation
  2. Plug this into the stress tensor equation and display the result in Cartesian
    and spherical coordinates?
Megan Ly
Is there a way to get Sympy to pull a constant inside of an integral without evaluating it? For example, rewriting c*Integral(f(x), (x,a,b)) as Integral(c*f(x), (x,a,b)). I'm struggling to get replace to work.
4 replies
Hello everyone!
Am Ruturaj, computer science student.
I want to contribute to this open source, so can anyone guide me please?
Cherish Sachdeva
Hello Everyone!
Hope you all are having a nice day.I want to integrate(cos(3theta)/(5-4cos(theta)),(theta,0,2pi)) this expression whose
answer is pi/12 but sympy does not executes the integral. So,I replaced cos(3
theta) with
4cos(theta)*3-3cos(theta) then I am obtaining -21pi/8.Also when I replaced cos(3theta) with
cos(theta) sympy executed the integral and gave -pi/2.Can anyone please elaborate why sympy is not executing
cos(3*theta) in integral?
Kalevi Suominen
Integrals of that type should be handled by the risch integrator but that part of the algorithm has not been implemented (yet). Currently, heurisch will be used as a fall-back integrator but it looks like there is a bug that leads to a wrong result.
Megan Ly
Is there a way to get the Latex printer to print Integral(x**2,(x,1,2)) as the string '\\int_{1}^{2} x^{2}\\, dx' instead of '\\int\\limits_{1}^{2} x^{2}\\, dx'?
Aaron Meurer
@meganly apparently latex(mode='inline') does that.
I wonder what the point of \limits is.
12x2dx\int_1^2 x^2\, dx vs. 12x2dx\int\limits_1^2 x^2\, dx
I guess \limits puts them directly above or below the integral sign. That does seem better for non-inline mode, especially if the limits can be long.
then again the non-limits version does look better even for equation mode when the limits are small. Maybe there should be some heuristics here.
Aaron Meurer
I opened sympy/sympy#21851 for this.
Megan Ly
Thanks. My reason for not wanting \limits is that the \limits command is not supported by MathQuill which causes issues for rendering Latex in a web browser.
hey hey, I saw in the release notes that ArraySymbol was added, however, I couldn't find the documentation on it, does anyone have any hints where to find more information?
Aaron Meurer
It looks like it isn't documented very much yet. You may have to take a look at the PR that added it for now to get some more information.
Thanks, I looked into the source and it basically is a way to conveniently get symbols with labes that make sense for arrays and array elements, but there is not other functionality that really requires documentation
Aaron Meurer
I think it should be roughly the same idea as MatrixSymbol
Any docs to build sympy from source?
2 replies
Lambdify uses _EvaluatorPrinter and I want it to behave slightly differently. Is there any option other than a lot of code duplication/modifying the source? Specifically, I don't want it to use return, but instead I want to use ASTs return (so that I can e.g. go for multiple lines in one go or use control flow logic from AST)
In particular, my solution would be to replace '''funcbody.append('return ({})'.format(self._exprrepr(expr)))''' with '''funcbody+=self._exprrepr(expr).splitlines()''' if expr is a CodeBlock on line 1092 of utilities/lambdify.py
also, gratz on the 1.9 release :D
Heyy devs, i am new to open source and i would certainly love to help in whichever way possible. Is there anything i could help you guys with?
When I want to copy some file from sympy and modify it slightly for my own library, can I do that? I would guess I need to add the Sympy license to my own project? If so, do I need to add all the "external" licenses of the sympy project too, or only the main sympy license?
Jason Ross
Is there a generalized Kronecker Delta in sympy?
@ThePauliPrinciple thanks for answering sir .
There is a levi-civita symbol, which you may be able to use @Jason-S-Ross
alternatively, you could explicitly write the pxp determinant https://docs.sympy.org/latest/tutorial/matrices.html#determinant
Jason Ross
@ThePauliPrinciple I think the levi-civita symbol would work. Good thinking
Hi Folks! Is it possible to use sympy to obtain the probability density function of a transformed random variable?
I am not expecting the results given by https://tinyurl.com/y2vzw6dp
Matías Senger

I am working in a simple QFT calculation and would like to do it using Sympy to learn (and also check my result). I have found the Quantum Mechanics module but cannot see how to start using it for my purpose. I have defined these quantities:

import sympy.physics.quantum as Q
vacuum = Q.OrthogonalKet(0)
annihilation_op = Q.Operator('a')
creation_op = Q.Dagger(annihilation_op)

and now I want to tell Sympy that an=nn1a |n\rangle = \sqrt{n} | n-1 \rangle if n>0n>0 else 00 and an=n+1n+1a^\dagger | n \rangle = \sqrt{n+1} | n+1 \rangle. How would I do this? Also, how do I impose the commutation relations between aa and aa^\dagger?


Hi all,

I am just a new user for SymPy. I am self learning this library for my undergrauate research. But in the middle of the process I am stucked with one code.

So I have defined a function with a subscript.

                     U_n= x^n + 1/x^n

When I consider (U_1)^3 I get (substitute n=1)

                    (U_1)^3 = (x+1/x)^3

Then after simplifying this I get

                    (U_1)^3 = (x^3 + 1/x^3) + 3(x+ 1/x)      

But one can see this answer as

                                   (U_1)^3 =   U_3 + 3U_1

How to get the output in terms of U_n 's ?

Can someone please give an idea how to build this code using SymPy. It would be a very big help for my research.

Thank you very much.

Gayanath Chandrasena.

You can use .subs to substitute x+1/x with U_1
But, you need a better definition of what you want
since U_1^3 is already in terms of U_n's?
Cheers! I'm having an issue with an integration where sympy outputs a piecewise with the condition dr*(k-1) < dr*(k+1), which is obviously true, since dr = symbols('Delta',real=True,positive=True,nonzero=True) and k = symbols('k',integer=True,real=True), yet sympy doesn't seem to thing this is the case. Am I doing something wrong?

Also I get

In [54]: print(ask(dr*(k-1) < dr*(k+1)))


In [55]: print(ask(dr*k-dr < dr*k+dr))
1 reply
So sympy seems to think that it can't multiply out the brackets for some reason?
I'm on sympy 1.8, btw.
@ThePauliPrinciple If I simplify U_1^3 I will be getting the answers in X terms.But if you further inspect that term can be converted again to U_n terms(ex: U_1^3=U_3+3U_1). Simple substitution does not work here. There should be a method to back-substitute.
But you still have not defined what it is you actually want, you only gave a single example, do you want an arbitrary expression to be expressed in terms of a sum_i a_i*U_i? with in the above case a_3=1 and a_1=3?

My function is
U_n= x^n + 1/x^n.

As an example when I compute (U_1)^3 I get (x^3 + 1/x^3) + 3(x+ 1/x).---------(i)

And if I compute (U_2)^2 I get x^4 + 1/x^4 + 2---------(ii)

But since in (i), (x^3 + 1/x^3)= U_3 and 3(x+ 1/x)=3U_1 I want to get the answer U_3+3U_1.

In (ii) I want to get the answer U_4 +2 since x^4 + 1/x^4 = U_4.

@GayanathSamuditha Although I dont know how you are going to interface with this, but I have created something which can do what you explained.
aman@amanUBUNTU:~/Desktop$ python3 Un.py 
Enter Value of n and k as in (U_n)^k
n = 1
k = 3
3*U_1 + U_3
aman@amanUBUNTU:~/Desktop$ python3 Un.py 
Enter Value of n and k as in (U_n)^k
n = 1 
k = 2
U_2 + 2
aman@amanUBUNTU:~/Desktop$ python3 Un.py 
Enter Value of n and k as in (U_n)^k
n = 2
k = 10
45*U_12 + 10*U_16 + U_20 + 210*U_4 + 120*U_8 + 252