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.

Reply via email to