Ondrej, Thanks for all your answers and your quick reply. I ask for forgiveness of my notably beginner's code developing, thanks for helping out! I attached a modified f.py that clarifies what I intend to integrate, along with an example. All and each of your observations are greatly appreciated, I'll be applying the guidelines from today on.
Maybe you want to know that what I'm developing is a code to symbolically perform Bayesian Statistics Inference for an environmental problem. In the future, I could put together some MLE, Bayeasian and MCMC algorithms to be included in Sympy, that I used to numerically implement in Matlab. Would it be a nice add-on? Thanks again, Angelica. ________________________________________ From: [EMAIL PROTECTED] [EMAIL PROTECTED] On Behalf Of Ondrej Certik [EMAIL PROTECTED] Sent: Tuesday, November 11, 2008 5:39 AM To: sympy-patches@googlegroups.com Cc: Echavarria Gregory, Maria Angelica Subject: Re: [PATCH] Patch is an added function Gaussian Integral in Integrals.py, to solve the definite integral of exp(poly). Hi Angelica! thanks very much for the patch and the work you did. I have couple comments below: On Tue, Nov 11, 2008 at 12:20 AM, <[EMAIL PROTECTED]> wrote: > > # HG changeset patch > # User Angelica Echavarria-Gregory Add your email address to your ~/.hgrc, e.g. mine looks like this: [ui] username = Ondrej Certik <[EMAIL PROTECTED]> in windows I guess it is some option of TortoiseHg. > # Date 1226355651 18000 > # Node ID 80e7a79d74d398ac99323fd100c56063f2bbedce > # Parent d66dc8f1f1077acc5776284d9fbbcc4c9871fb01 > The function to PLEASE MERGE WITH ALL OF SYMPY FOR IT TO WORK ALONG WITH > OTHER INTEGRALS; > # Please, let's implement this function, we need Sympy to, for example, > recognize when he find this expresion next to a different expression and be > able to solve all together: > > # Thanks Stephan who came up with the idea of converting the poly into a > Gaussian Integral. > > diff -r d66dc8f1f107 -r 80e7a79d74d3 sympy/integrals/integrals.py > --- a/sympy/integrals/integrals.py Fri Oct 24 19:09:04 2008 +0200 > +++ b/sympy/integrals/integrals.py Mon Nov 10 17:20:51 2008 -0500 > @@ -376,4 +376,28 @@ > if isinstance(integral, Integral): > return integral.doit() > else: > - return integral > + return integral > + > + > +# The function to PLEASE MERGE WITH ALL OF SYMPY FOR IT TO WORK ALONG WITH > OTHER INTEGRALS; > +# Please, let's implement this function, we need Sympy to, for example, > recognize when he find this expresion next to a different expression and be > able to solve all together: > + > +# Thanks Stephan who came up with the idea of converting the poly into a > Gaussian Integral. > + > +import sympy > +from sympy import Symbol, integrate, sqrt, simplify, exp, erf, pprint If you use the import statements like this, it won't import: $ python Python 2.5.2 (r252:60911, Sep 29 2008, 21:15:13) [GCC 4.3.2] on linux2 Type "help", "copyright", "credits" or "license" for more information. >>> import sympy Traceback (most recent call last): File "<stdin>", line 1, in <module> File "sympy/__init__.py", line 28, in <module> from integrals import * File "sympy/integrals/__init__.py", line 8, in <module> from integrals import integrate, Integral File "sympy/integrals/integrals.py", line 388, in <module> from sympy import Symbol, integrate, sqrt, simplify, exp, erf, pprint ImportError: cannot import name integrate Also, the convention is to put all import statements on top of the module. > +def GaussIntegral(expr,param,lowerbound,upperbound,z=Symbol('z')): We use Python Style Guidelines: http://www.python.org/dev/peps/pep-0008/ so the function should be namedg auss_integral() and also please use spaces after commans, e.g.: gauss_integral(expr, param, lowerbound, upperbound, z=Symbol('z')) > + """Finds the definite integral of exp(second order polynomial) > expression trough transforming into Gaussian Integral""" Could you please also put some example in the docstring how it is supposed to be used? > + exponent = expr.args[0] > + coeff = exponent.as_poly(param).coeffs > + t1 = sqrt(coeff[0])*param > + t2 = coeff[1]*param/(2*t1) > + square = t1+t2 > + remainder = simplify(exponent-square**2) > + jacobn = 1/square.diff(param) > + result = integrate(exp(-z**2),z).subs(z,square)*jacobn*exp(remainder) > + result The last line (result) is not necessary, right? > + > A=(result.subs(param,upperbound).evalf())-(result.subs(param,lowerbound).evalf()) > + return pprint(A) pprint returns None -- it prints directly to the terminal. I suggest you simply return A, and handle the pprint outside of this function, e.g. the user himself. I am trying to run it, but I am getting some problems -- if you could please give me some examples how you run the function, it be cool. I tried: x = Symbol("x") gauss_integral(exp(x**2), x, -oo, oo) and this is what I get: File "f.py", line 20, in <module> gauss_integral(exp(x**2), x, -oo, oo) File "f.py", line 9, in gauss_integral t2 = coeff[1]*param/(2*t1) IndexError: tuple index out of range See the attached f.py. If you could please modify the f.py to run together with some examples, I'll then help you integrate it well with sympy. Thanks again for your work! Ondrej --~--~---------~--~----~------------~-------~--~----~ You received this message because you are subscribed to the Google Groups "sympy-patches" group. To post to this group, send email to sympy-patches@googlegroups.com To unsubscribe from this group, send email to [EMAIL PROTECTED] For more options, visit this group at http://groups.google.com/group/sympy-patches?hl=en -~----------~----~----~----~------~----~------~--~---
#This function solves definite integrals subject to: # 1. 'expr' is of the form ((x-a)**2-(y-b)**2)/(c*Z), where x, y and z are symbols and a, b, c constants. Any constants can be added. # 2. integrals such as integral(exp(x**2)) are already solved by sympy, so this is not a good example. # 3. 'lowerbound' and 'upperbound' are reals different than -oo and oo (BUT I WOULD LIKE and need THE FUNCTION TO HANDLE -OO AND OO, IT'S JUST THAT I DON'T KNOW HOW TO MAKE THAT POSSIBLE) # 4. 'param' could be either of the symbols in the 'expr'. (Actually I need to integrate with respect to all of them, in a sequence) # 5. Capability to add: the 'result' should be again integrated with respect to the following symbol, and so on, like when using integrate(integrate(integrate(f(x,y,z)x).doit(),y).doit(),z).doit(), with definite bounds for each parameter. # 6. Please see example of how I will use it below the function. import sympy from sympy import Symbol, integrate, sqrt, simplify, exp, erf, pprint, oo def gauss_integral(expr, param, lowerbound, upperbound, z=Symbol('z')): """Finds the definite integral of exp(second order polynomial) expression trough transforming into Gaussian Integral""" exponent = expr.args[0] coeff = exponent.as_poly(param).coeffs t1 = sqrt(coeff[0])*param t2 = coeff[1]*param/(2*t1) square = t1+t2 remainder = simplify(exponent-square**2) jacobn = 1/square.diff(param) result = integrate(exp(-z**2),z).subs(z,square)*jacobn*exp(remainder) A=(result.subs(param,upperbound).evalf())-(result.subs(param,lowerbound).evalf()) return A #x = Symbol("x") #gauss_integral(exp(x**2), x, -oo, oo) # This is not required but sure useful, but it doesn't work for me either with -oo and oo. # This is already implemented in Sympy when the bonds are not oo or -oo, I get this: >>> integrate(exp(x**2),x).doit() -I*pi**(1/2)*erf(I*x)/2 # EXAMPLE --THIS IS PART OF THE CODE I'M TRYING TO DEVELOP-- This is what I get when I run it. # Indexes can take any value within range of arrays in the real copy, but were replaced by cnstants in the example. Please disregard unnecessary terms. >>> execfile('/Python25/Scripts/Gaussian Integral.py') # Sure I'll use the >>> guidelines freom now, thanks for the observation!! >>> import numpy >>> from numpy import array, average, std, cos, sin, pi >>> from __future__ import division >>> from sympy import Symbol, sqrt, simplify, integrate, exp, erf, pprint, pi >>> t=numpy.array([1.0, 1.1, 1.25, 1.5, 3.0, 5.0, 6.0, 8.0, 11.0]) >>> x=numpy.array([1.00,2.00,3.00,4.00,5.00,6.00,7.00,8.00,9.00,10.00,11.00,12.00,13.00,14.00,15.00]) >>> y=numpy.array([3.00,4.00,5.00,6.00,7.00,8.00,9.00,10.00,11.00,12.00,13.00,14.00,15.00,16.00,17.00]) >>> x0=0. >>> y0=0. >>> vx=Symbol('vx') >>> vy=Symbol('vy') >>> D=Symbol('D') >>> expr = >>> exp(((-(x[0]-x0-(vx*t[0])))**2)+(-(y[0]-y0-(vy*t[0]))**2)/(4.00*D*t[0])) >>> GaussIntegral(expr,vx,1,10) # Any boundaries, this is just an example. >>> Wonderful if it worked with infinites: like gauss_integral(expr,vx,0, oo). 2 2.25 1.5*vy 0.25*vy - ---- + ------ - -------- D D D 0.886226925452758*e >>> # Nice result! # I will not pprint from code for now on. I am learning to be a good pupil. Just to show you what I get. # Following, I'd need to partially integrate this result with respect to vy, and the following answer, with respect to D; but LOOK! the expression is no longer of the 'expr' type, # it is multiplied by a constant, so I need to make use of other functions of Sympy to recognize this and do what necessary to separete # these two factors and get the correct integration.