Search code examples
pythonscipynumerical-integration

Integrate in a if statement? (python)


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)

Solution

  • 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.)