Symbolic computation with partial derivatives and a double integral

  • Thread starter Thread starter fluidistic
  • Start date Start date
  • Tags Tags
    Calculus Python
AI Thread Summary
The discussion revolves around computing an integral involving a vector field in polar coordinates within a quarter annulus. The user has derived the x-component of the vector field and applied the chain rule to express the derivative with respect to x. They implemented a Simpy code to compute the integral but are seeking further simplification of the resulting expression. Suggestions include substituting logarithmic terms and ensuring the variables are defined as real and positive to facilitate simplification. Comparisons with Mathematica's output reveal potential discrepancies and opportunities for further refinement of the expression.
fluidistic
Gold Member
Messages
3,928
Reaction score
272
I have an expression of a vector field ##\vec J = J_r \hat r + J_\theta \hat \theta## in polar coordinates, in a region that's a quarter of an annulus, where theta ranges between 0 and pi/2 and r goes from ri to ro.
I need to compute ##-\int _\Omega T(S_{xx}-S_{yy}) \frac{\partial J_x}{\partial x}d\Omega =-(S_{xx}-S_{yy})\int_0^{\pi/2} \int_{ri}^{ro} [T_0 + \Delta T \ln(r/ro) / \ln(ro/ri)] \frac{\partial J_x}{\partial x} rdrd\theta## where all the parameters in the last expressions are known constants.
The thing is that I do not have the ##\hat x## component of ##\vec J##, I have the ##\hat r## and ##\hat \theta## components only.
Therefore I need to obtain ##J_x## first. I get that ##J_x = J_r \cos(\theta)-J_\theta \sin(\theta)##.
Using the chain rule I get that ##\frac{\partial }{\partial x}=\cos(\theta)\frac{\partial}{\partial r}-\frac{\sin(\theta)}{r}\frac{\partial}{\partial \theta}##.
Having derived this, I can now write the Simpy code that should give me the sought integral.

Python:
import simpy as sp

# Define the variables
r, theta, H0, ri, ro, rho, T_0, Delta_T, S_xx, S_yy = sp.symbols('r theta H0 ri ro rho T_0 Delta_T S_xx S_yy')

# Define H0
H0_expr = -(S_xx - S_yy) / (4 * rho) * Delta_T / sp.log(ro/ri)

# Define the components of the vector field J_r and J_theta
J_r = 2 * H0 * (1/r - 1/(ri**2 + ro**2)*(r + ri**2*ro**2/r**3)) * sp.cos(2*theta)
J_theta = H0 * (1/(ri**2 + ro**2)*(2*r - 2*ri**2*ro**2/r**3)) * sp.sin(2*theta)

J_x = J_r * sp.cos(theta) - J_theta * sp.sin(theta)
dJ_x_dx = sp.diff(J_x, r) * sp.cos(theta) - 1/r * sp.sin(theta) * sp.diff(J_x, theta) * sp.sin(theta)

T = T_0 + Delta_T * sp.log(r/ro) / sp.log(ro/ri)

# Set up the integral for Q
Q = sp.integrate(r * T * (S_xx - S_yy) * dJ_x_dx, (theta, 0, sp.pi/2), (r, ri, ro))

# Substitute H0 in Q
Q_sub = Q.subs(H0, H0_expr)
Q_sub = sp.simplify(Q_sub)
Q_sub = sp.collect(Q_sub, [ri, ro, sp.log(ro/ri)])
Q_sub = sp.factor(Q_sub)
print(f"This is Q: {Q}")

The code yields
This is Q: -Delta_T*(S_xx - S_yy)**2*(-16*Delta_T*ri**2*log(ri/ro)**2 + 15*pi*Delta_T*ri**2*log(ri/ro)**2 + 80*Delta_T*ri**2*log(ri/ro) + 30*pi*Delta_T*ri**2*log(ri/ro) - 64*Delta_T*ri**2 - 15*pi*Delta_T*ri**2 - 16*Delta_T*ro**2*log(ri/ro)**2 + 15*pi*Delta_T*ro**2*log(ri/ro)**2 + 48*Delta_T*ro**2*log(ri/ro) + 15*pi*Delta_T*ro**2 + 64*Delta_T*ro**2 - 32*T_0*ri**2*log(ri)*log(ro/ri) + 30*pi*T_0*ri**2*log(ri)*log(ro/ri) - 30*pi*T_0*ri**2*log(ro)*log(ro/ri) + 32*T_0*ri**2*log(ro)*log(ro/ri) + 32*T_0*ri**2*log(ro/ri) + 30*pi*T_0*ri**2*log(ro/ri) - 32*T_0*ro**2*log(ri)*log(ro/ri) + 30*pi*T_0*ro**2*log(ri)*log(ro/ri) - 30*pi*T_0*ro**2*log(ro)*log(ro/ri) + 32*T_0*ro**2*log(ro)*log(ro/ri) - 30*pi*T_0*ro**2*log(ro/ri) - 32*T_0*ro**2*log(ro/ri))/(480*rho*(ri**2 + ro**2)*log(ro/ri)**2).

I wonder what Mathematica (or a better Simpy code) can yield, I strongly suspect the expression can be further simplified. Can someone confirm my result or show me a shoter one, which I could then compare?
 
Physics news on Phys.org
It is apparent from the final expression for Q that SymPy is not eliminating one of \log(r_o/r_i) and \log(r_i/r_o) in favour of the other, so you could try explicitly telling it to do that by making the substitution
Python:
Q_sub = Q_sub.subs(sp.log(ri/ro), -sp.log(ro/ri))
before simplifying and collecting terms.

If you explicitly tell SymPy that r_o and r_i are real and positive by declaring them as
Python:
ro = sp.Dummy('ro', real=True, positive=True)
etc. rather than just using sp.symbols, then SymPy might be able to apply further simplifications which don't hold for negative or complex arguments.
 
Last edited:
  • Informative
  • Like
Likes fluidistic and renormalize
fluidistic said:
I wonder what Mathematica (or a better Simpy code) can yield, I strongly suspect the expression can be further simplified. Can someone confirm my result or show me a shoter one, which I could then compare?
Here is what Mathematica gives for your heat expression:
1733175460483.png

Cleaning up this result a bit, I repeat it below in LaTeX, along with the Joule heat from your previous thread https://www.physicsforums.com/threads/evaluating-integrals-derivatives-etc-with-simpy.1067101/ :$$Q_{\text{Joule}}=\frac{\pi\,\left(\Delta T\right)^{2}\left(S_{xx}-S_{yy}\right)^2\left[\left(r_{i}^{2}+r_{o}^{2}\right)\log\left(\frac{r_{o}}{r_{i}}\right)+r_{i}^{2}-r_{o}^{2}\right]}{16\,\rho\left(r_{i}^{2}+r_{o}^{2}\right)\log^{2}\left(\frac{r_{o}}{r_{i}}\right)}$$$$Q_{\text{Bridgman}}=\frac{\pi\,\Delta T\left(S_{xx}-S_{yy}\right)^2\left[\Delta T\left(r_{o}^{2}-r_{i}^{2}\right)-2\log\left(\frac{r_{o}}{r_{i}}\right)\left(\Delta T\,r_{i}^{2}+T_{0}\left(r_{o}^{2}-r_{i}^{2}\right)\right)\right]}{16\,\rho\left(r_{i}^{2}+r_{o}^{2}\right)\log^{2}\left(\frac{r_{o}}{r_{i}}\right)}$$Perhaps you can use this to compare to your SymPy expression after applying the suggestions of @pasmith in post #2?
 
Last edited:
pasmith said:
It is apparent from the final expression for Q that SymPy is not eliminating one of \log(r_o/r_i) and \log(r_i/r_o) in favour of the other, so you could try explicitly telling it to do that by making the substitution
Python:
Q_sub = Q_sub.subs(sp.log(ri/ro), -sp.log(ro/ri))
before simplifying and collecting terms.

If you explicitly tell SymPy that r_o and r_i are real and positive by declaring them as
Python:
ro = sp.Dummy('ro', real=True, positive=True)
etc. rather than just using sp.symbols, then SymPy might be able to apply further simplifications which don't hold for negative or complex arguments.
Thanks for the suggestion. The expression is still very large after the substitution of the logs.
And the computations is taking forever (I will have to abort) when I specify both ri and ro to be real and positive.
 
renormalize said:
Here is what Mathematica gives for your heat expression:
View attachment 354047
Cleaning up this result a bit, I repeat it below in LaTeX, along with the Joule heat from your previous thread https://www.physicsforums.com/threads/evaluating-integrals-derivatives-etc-with-simpy.1067101/ :$$Q_{\text{Joule}}=\frac{\pi\,\left(\Delta T\right)^{2}\left(S_{xx}-S_{yy}\right)\left[\left(r_{i}^{2}+r_{o}^{2}\right)\log\left(\frac{r_{o}}{r_{i}}\right)+r_{i}^{2}-r_{o}^{2}\right]}{16\,\rho\left(r_{i}^{2}+r_{o}^{2}\right)\log^{2}\left(\frac{r_{o}}{r_{i}}\right)}$$$$Q_{\text{Bridgman}}=\frac{\pi\,\Delta T\left(S_{xx}-S_{yy}\right)\left[\Delta T\left(r_{o}^{2}-r_{i}^{2}\right)-2\log\left(\frac{r_{o}}{r_{i}}\right)\left(\Delta T\,r_{i}^{2}+T_{0}\left(r_{o}^{2}-r_{i}^{2}\right)\right)\right]}{16\,\rho\left(r_{i}^{2}+r_{o}^{2}\right)\log^{2}\left(\frac{r_{o}}{r_{i}}\right)}$$Perhaps you can use this to compare to your SymPy expression after applying the suggestions of @pasmith in post #2?
That's amazing, thanks a lot! Woah... There's a little typo in your latex (missing square for the S_xx and S_yy term).
This helps a lot.
 
  • Like
Likes renormalize
fluidistic said:
There's a little typo in your latex (missing square for the S_xx and S_yy term).
Thanks, fixed it in post #3.
 

Similar threads

Back
Top