Hi everyone, I am new to Sage and coding, and I am trying to solve a tricky integral from equation 13 in this physics paper: www.nature.com/articles/srep11876
Here is the integrand in code form (I would upload a screenshot of the full equation, but I do not have enough karma points):
sage: var ('alpha beta E_o a d_o z_o n u')
(alpha, beta, E_o, a, d_o, z_o, n, u)
sage: j = CC.0
sage: integrand = e**(-j * n * u) / (((a + d_o + 0.5 * z_o * (((16 * pi) / (alpha * beta)) ** (-3))) + 0.5 * z_o * (((16 *
....: pi) / (alpha * beta)) ** (-3)) * sin(u)) ** 3) - 1
The integral is part of a Fourier coefficient, so it is a definite integral with limits -pi and pi. I want to integrate with respect to u and then solve the full equation for beta. All of the other variables are known.
However, when I try to integrate, I get an error:
sage: integrate(integrand, u, -pi, pi)
---------------------------------------------------------------------------
RuntimeError Traceback (most recent call last)
<ipython-input-32-84f9192f9549> in <module>()
----> 1 integrate(integrand, u, -pi, pi)
/opt/sagemath-9.1/local/lib/python3.7/site-packages/sage/misc/functional.py in integral(x, *args, **kwds)
751 """
752 if hasattr(x, 'integral'):
--> 753 return x.integral(*args, **kwds)
754 else:
755 from sage.symbolic.ring import SR
/opt/sagemath-9.1/local/lib/python3.7/site-packages/sage/symbolic/expression.pyx in sage.symbolic.expression.Expression.integral (build/cythonized/sage/symbolic/expression.cpp:64474)()
12434 R = ring.SR
12435 return R(integral(f, v, a, b, **kwds))
> 12436 return integral(self, *args, **kwds)
12437
12438 integrate = integral
/opt/sagemath-9.1/local/lib/python3.7/site-packages/sage/symbolic/integration/integral.py in integrate(expression, v, a, b, algorithm, hold)
929 return indefinite_integral(expression, v, hold=hold)
930 else:
--> 931 return definite_integral(expression, v, a, b, hold=hold)
932
933
/opt/sagemath-9.1/local/lib/python3.7/site-packages/sage/symbolic/function.pyx in sage.symbolic.function.BuiltinFunction.__call__ (build/cythonized/sage/symbolic/function.cpp:12262)()
1136 res = self._evalf_try_(*args)
1137 if res is None:
-> 1138 res = super(BuiltinFunction, self).__call__(
1139 *args, coerce=coerce, hold=hold)
1140
/opt/sagemath-9.1/local/lib/python3.7/site-packages/sage/symbolic/function.pyx in sage.symbolic.function.Function.__call__ (build/cythonized/sage/symbolic/function.cpp:6938)()
595 for i from 0 <= i < len(args):
596 vec.push_back((<Expression>args[i])._gobj)
--> 597 res = g_function_evalv(self._serial, vec, hold)
598 elif self._nargs == 1:
599 res = g_function_eval1(self._serial,
/opt/sagemath-9.1/local/lib/python3.7/site-packages/sage/symbolic/function.pyx in sage.symbolic.function.BuiltinFunction._evalf_or_eval_ (build/cythonized/sage/symbolic/function.cpp:13412)()
1224 res = self._evalf_try_(*args)
1225 if res is None:
-> 1226 return self._eval0_(*args)
1227 else:
1228 return res
/opt/sagemath-9.1/local/lib/python3.7/site-packages/sage/symbolic/integration/integral.py in _eval_(self, f, x, a, b)
203 for integrator in self.integrators:
204 try:
--> 205 A = integrator(*args)
206 except (NotImplementedError, TypeError, AttributeError):
207 pass
/opt/sagemath-9.1/local/lib/python3.7/site-packages/sage/symbolic/integration/external.py in maxima_integrator(expression, v, a, b)
44 result = maxima.sr_integral(expression,v)
45 else:
---> 46 result = maxima.sr_integral(expression, v, a, b)
47 return result._sage_()
48
/opt/sagemath-9.1/local/lib/python3.7/site-packages/sage/interfaces/maxima_lib.py in sr_integral(self, *args)
788 """
789 try:
--> 790 return max_to_sr(maxima_eval(([max_integrate],[sr_to_max(SR(a)) for a in args])))
791 except RuntimeError as error:
792 s = str(error)
/opt/sagemath-9.1/local/lib/python3.7/site-packages/sage/libs/ecl.pyx in sage.libs.ecl.EclObject.__call__ (build/cythonized/sage/libs/ecl.c:7794)()
803 """
804 lispargs = EclObject(list(args))
--> 805 return ecl_wrap(ecl_safe_apply(self.obj,(<EclObject>lispargs).obj))
806
807 def __richcmp__(left, right, int op):
/opt/sagemath-9.1/local/lib/python3.7/site-packages/sage/libs/ecl.pyx in sage.libs.ecl.ecl_safe_apply (build/cythonized/sage/libs/ecl.c:5456)()
375 if ecl_nvalues > 1:
376 s = si_coerce_to_base_string(ecl_values(1))
--> 377 raise RuntimeError("ECL says: {}".format(
378 char_to_str(ecl_base_string_pointer_safe(s))))
379 else:
RuntimeError: ECL says: 1.8189894035458565e-12 is not of type (UNSIGNED-BYTE 15).
I have no idea where 1.8189894035458565e-12 came from, and I get the same error after assigning values to all of the variables except beta and u.
Has anyone encountered this error before? Is there a better way I should evaluate this?
