Minimize constant of a definite integral using Scipy

Viewed 78

I am trying to optimize the following problem using Scipy

formula

Currently I have defined my equations and attempted to use the quad function and minimize to solve the equation.

import numpy as np
from scipy.integrate import quad
from scipy.optimize import minimize


def g_x(x, cv):
    cv2 = np.square(cv)
    return (x+1-np.sqrt(np.square(x)+2*cv2*x+1))/(1-cv2)

def g_approx(x, c):
    return (c*x)/(1+c*x)

def integrand(x, cv, c):
    return g_x(x, cv)-g_approx(x, c)

def minimization_function(c, cv, xm):
    return quad(integrand, 0, xm, args=c)

c_0 = 1
cv = 0.6
xm = 15
res = minimize(minimization_function, c_0, method='CG', tol=1.e-2, options={'gtol': 0.01, 'maxiter': 5})

When trying to solve this problem I get the following error.

scipy\optimize\optimize.py:261, in _prepare_scalar_function(fun, x0, jac, args, bounds, epsilon, finite_diff_rel_step, hess)
257     bounds = (-np.inf, np.inf)
259 # ScalarFunction caches. Reuse of fun(x) during grad
...
132 # Overwriting results in undefined behaviour because
133 # fun(self.x) will change self.x, with the two no longer linked.
--> 134 return fun(np.copy(x), *args)

TypeError: minimization_function() missing 2 required positional arguments: 'cv' and 'xm'

So I know at some point I need to give C_v and X_m as they are both missing to solve the problem, however, I'm unsure of how to do this. Could anyone help me understand how to define variables and constants in these equations? I've tried looking through previous answers, but I can't seem to make them work.

1 Answers

Both quad and minimize expect a function with signature fun(x, *args). Thus, you either need to wrap the function by a lambda expression or use the args parameter. Note also that quad doesn't return a scalar, so you need to extract the integral value from the returned tuple:

# ... your other functions and constants here ...

def minimization_function(c, cv, xm):
    return quad(lambda x: integrand(x, cv, c), 0, xm)[0]

res = minimize(lambda c: minimization_function(c, cv, xm), c_0)

Using the args parameter instead, your code would look like this:

# ... your other functions and constants here ...

def minimization_function(c, cv, xm):
    return quad(integrandx, 0, xm, args=(cv, c))[0]

res = minimize(minimization_function, c_0, args=(cv, c))

PS: There's no need to evaluate the integral numerically by quad. Instead, you can simply use the closed-form expression (see WolframAlpha if you don't want to derive it on your own).

Related