pythonfor-loopsympybessel-functions

How to Define a Piecewise Function in Python with Multiple Variables


I am trying to develop a plot for my helioseismology class and the question had provided a piecewise function describing the dynamics of the "fluids" in a star as if it is one thing its this and if its another its that. I am receiving over and over again this 'Mul' object cannot be interpreted as an integer but I am working with numbers in the reals not just the integer set. I do not know how to get around this and need guidance. The code is as follows.

import sympy as sy
from sympy import *
from sympy.physics.units import Unit
import numpy as np
import sys
import math
import scipy as sp
from scipy import special

phi = Symbol('phi', Variable = True)
x = Symbol('x', Variable = True, Real = True)
t = Symbol('t', Variable = True, Real = True)
xi = Symbol('xi', Function = True)
Solar_Radius = Symbol('R', Constant = True, unit = "meters")
Sound_Speed = Symbol('c', Constant = True, unit = "meters per second", Real = True)
gamma = Symbol('gamma', Constant = True)
gravity = Symbol('g', Constant = True, unit = "meters per second per second")

Solar_Radius = 6.963 * 10 ** 6
gamma = 5/3
g = 274.8265625336
gas_constant = 8201.25
c = 8.1 * 10 ** 3

for t in range(0,x/c):
    xi[x,t] = 0
for t in range(x/c,00):
    xi[x,t] = (1/2)*sy.exp(gamma*g*x/(2*c**2))*mpmath.besselj(0, (gamma*g/(2*c)*sy.sqrt(t**2 - ((x/c)**2))),derivative = 0)

Full Traceback:

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-50-3506376f1686> in <module>()
----> 1 for t in range(0,x/c):
      2     xi[x,t] = 0
      3 for t in range(x/c,00):
      4     xi[x,t] = (1/2)*sy.exp(gamma*g*x/(2*c**2))*mpmath.besselj(0, (gamma*g/(2*c)*sy.sqrt(t**2 - ((x/c)**2))),derivative = 0)

TypeError: 'Mul' object cannot be interpreted as an integer

Solution

  • There are quite a few issues here:

    What you want here is Piecewise. The syntax is Piecewise((expr, cond), (expr, cond), ..., (expr, True)), where expr is an expression that is used when cond is true ((expr, True) is the "otherwise" condition).

    For your example, I believe you want

    expr = Piecewise((0, t < x/c), (sy.exp(gamma*g*x/(2*c**2))*sy.besselj(0, (gamma*g/(2*c)*sy.sqrt(t**2 - (x/c)**2)))/2, t >= x/c))
    

    If you want to turn this into a numeric function in x and t, use

    xi = lambdify((x, t), expr)