[Python] Is this solvable using sympy? - no output after 1 day

The original question:

My approach: The electric field is already given in the textbook as

\begin{equation}
E(x, y, z) = \frac{q}{4\pi\epsilon_0} \left[
\frac{\hat{x}x + \hat{y}y + \hat{z}(z - z_0)}{(x^2 + y^2 + (z - z_0)^2)^{3/2}} +
\frac{\hat{x}x + \hat{y}y + \hat{z}(z + z_0)}{(x^2 + y^2 + (z + z_0)^2)^{3/2}}
\right]
\end{equation}


So I used this equation in the textbook:

\begin{equation}
U=\frac{\epsilon_0}{2}\int E^2d^3\!x=
\end{equation}


To calculate this, I asked GPT-4 to write a python program to compute this using sympy. Below is the code:

from sympy import *

# Define the symbols
x, y, z = symbols('x y z')
q = symbols('q', real=True, positive=True)
epsilon_0 = symbols('epsilon_0', real=True, positive=True)
z0 = symbols('z0', real=True, positive=True)

# Electric field components due to the image charges
E1_x = q / (4 * pi * epsilon_0) * x / ((x**2 + y**2 + (z - z0)**2)**(3/2))
E1_y = q / (4 * pi * epsilon_0) * y / ((x**2 + y**2 + (z - z0)**2)**(3/2))
E1_z = q / (4 * pi * epsilon_0) * (z - z0) / ((x**2 + y**2 + (z - z0)**2)**(3/2))

E2_x = q / (4 * pi * epsilon_0) * x / ((x**2 + y**2 + (z + z0)**2)**(3/2))
E2_y = q / (4 * pi * epsilon_0) * y / ((x**2 + y**2 + (z + z0)**2)**(3/2))
E2_z = -q / (4 * pi * epsilon_0) * (z + z0) / ((x**2 + y**2 + (z + z0)**2)**(3/2))

# Total electric field components
E_x = E1_x + E2_x
E_y = E1_y + E2_y
E_z = E1_z + E2_z

# Squaring the electric field components and summing them
E_squared = E_x**2 + E_y**2 + E_z**2

# Integrate E_squared over the half-space where z > 0
U = (epsilon_0 / 2) * integrate(integrate(integrate(E_squared, (x, -oo, oo)), (y, -oo, oo)), (z, 0, oo))

U.simplify()


I ran the code and expected to get a result in some day. I've left it running for around 1-2 days and there is no output. It is constantly using 100% CPU (I have 8P+2E cores). I'm not sure if this is not solvable in sympy or is it just really complicated.

Continue reading...