0

Is it possible for me to integrate but set my integration limits in a loop? So, for instance, integrate my function from 8 to 8.9 until it reaches a value of 2.5? Thank you!

f1 = lambda x,c : c/x
integrate1=quad(f, 8,8.97, args=c1)
print(integrate1)

Maybe?

for index in range(8,9):
    f1 = lambda x,c: c/x
    integrate1 = quad(f, index, index+0.1, args=c1)
    print(integrate1)
jayelm
  • 7,236
  • 5
  • 43
  • 61
  • 1
    http://docs.python.org/release/2.5.1/ref/indentation.html – Sinkingpoint Jan 20 '14 at 21:04
  • In the lower code block are several syntax errors, after `for ...` there should be a `:` and the following inner part of the loop should be indented. Also your name `f` and `f1` should probably be the same. –  Jan 20 '14 at 21:04
  • I'm a little confused on what you're asking here, what do you mean by " integrate my function from 8 to 8.9 until it reaches a value of 2.5"? – wnnmaw Jan 20 '14 at 21:15
  • Well, I have a distribution c/x and essentially I need to find the 95th percentile (confidence interval for that distribution). The lower end of my iteration was 8 and the higher, 20. So, I need to integrate from, say to to 8.5 to get a confidence interval of 2.5. I was trying to somehow loop the integral so that when it found that the sum was 2.5, it would stop. Hopefully, that's a bit clearer! – Arlena Sofia Jan 20 '14 at 22:43
  • Maybe I can also do a sum instead...say integral of c/x = c*ln(x). I will try to do it this way as well as see what happens... say, for index in range(small steps): res = res c*ln(index) if res <= 2.5: res = res – Arlena Sofia Jan 20 '14 at 22:44

2 Answers2

1

for a non-integer loop with a fixed step you may do something like that:

for val in xrange(80, 90):
    val /= 10
    ........

or

val = 8.0
while val<8.9:
    <do your worst>
    val += step

Rounding on the way may be good idea - to get required precision

volcano
  • 3,578
  • 21
  • 28
0

Well, obviously you can do that:

import numpy as np
import scipy.integrate as si

def test_fn(x, c):
    return c / x

def main():
    lower_limit =  8.0
    target_area =  2.5
    c_val       = 42.0
    for upper_limit in np.arange(8., 10., 0.1):
        area = si.quad(test_fn, lower_limit, upper_limit, args=(c_val,))
        if area >= target_area:
            print("Area is {} at ul={}".format(area, upper_limit))
            break

if __name__=="__main__":
    main()

but your step interval limits your result accuracy, and you're doing an awful lot of unnecessary calculations (== slow).

As @Jakob_Weisblat says, you can switch to a binary search algorithm. That's faster, but you have to do some extra bookkeeping. Why not delegate?

I suggest turning it into a metaproblem: solve for the upper limit such that integrating results in your desired target value:

import functools
import scipy.integrate as si
import scipy.optimize as so

def test_fn(x, c):
    return c / x

def integrate_fn(ul, fn, ll, args, target):
    return si.quad(fn, ll, ul, args=args) - target

def main():
    lower_limit =  8.0
    target_area =  2.5
    c_val       = 42.0
    sol = so.brentq(
        integrate_fn, lower_limit, 20.0, args=(
            test_fn, lower_limit, (c_val,), target_area
        )
    )
    print(sol)

if __name__=="__main__":
    main()

(Do note that this code is untested, as this machine does not have scipy installed.)

Hugh Bothwell
  • 55,315
  • 8
  • 84
  • 99