Find limits satisfying a condition on an integral

I have some data that I have interpolated with even spacing from measured data. I can integrate this between defined limits etc using trapz. However, I want to turn the problem on its head and find for what limit does the result of definite integral satisfy a condition. As an example:

x = linspace(0, 2*pi, 20);
y = sin(x);

lowLimIndex = 1;
upLimIndex = find(x>pi,1);
halfInt = trapz(x(lowLimIndex:upLimIndex),y(lowLimIndex:upLimIndex))

halfInt = 1.9682

Or close enough to the known result of 2... But in my case, I effectively want to find the value of x that results in the integral = 2.

I could fit an equation to my data and use solve to solve it in symbolic mode I suppose.... Seems there ought to be a way to do this using the optimization tools though?

Any help appreciated!

[big ted]

Thanks all. Granted, using a periodic function in my example might not have been the best choice! To be clear, my actual data isn't periodic, but this gives me a good footing to start from. For some reason, I've never quite gotten my head around the use of function handles as Torsten uses in his proof. This is a good reminder for me to knuckle down and do so!

Answer 1

[Torsten]

x = linspace(0, 2*pi, 20);
y = sin(x);
lowLimIndex = 1;
target = 2;
[~,upLimIndex] = min(abs(cumtrapz(x(lowLimIndex:end),y(lowLimIndex:end))-target));
x(upLimIndex)

ans = 2.9762

integral(@(x)sin(x),x(lowLimIndex),x(upLimIndex))

ans = 1.9864

Answer 2

[John D'Errico]

You have posed a interesting problem, with a subtle flaw that possibly you don't understand.

You asked to find an integral limit such that integral(sin(x)) == 2.

Ok. first, analytically, this might seem pretty easy.

syms x

I'll pick some arbitrary lower limit, not random, since I want it to be repeatable.

lolim = 0.5;
syms uplim

F = sin(x);

upsol = solve(int(F,[lolim,uplim]) == 2,returnconditions = true)

upsol = struct with fields:
         uplim: [2x1 sym]
    parameters: k
    conditions: [2x1 sym]

disp(char(upsol.uplim))

[acos(cos(1/2) - 2) + 2*k*pi; 2*k*pi - acos(cos(1/2) - 2)]

So there are infinitely many solutions. A fundamental one happens when k == 0.

disp(char(subs(upsol.uplim,'k',0)))

[acos(cos(1/2) - 2); -acos(cos(1/2) - 2)]

double(subs(upsol.uplim,'k',0))

ans = 
   3.1416 - 0.4899i
  -3.1416 + 0.4899i

Do you see? The solution is complex. Effectively, there are no solutions to the problem in real numbers.

If you set the lower limit to 0, then things get better.

int(sin(x),[0,pi])

ans = 

2

The problem is, if you set the lower limit to something other than some integer multiple of 2*pi, then you can't achieve an integral of 2.

And of course, it gets worse yet, if you use trapz. The problem is, no matter how small you st the increment, trapz will ALWAYS under-estimate the integral of a purely convex function like sin(x). For example...

fun = @sin;
X = linspace(0,pi,10);
format long g
trapz(X,fun(X))

ans = 
          1.97965081121648

X = linspace(0,pi,100);
trapz(X,fun(X))

ans = 
          1.99983216389399

What happens ism no matter how small the stepsize, trapz will underestimate the integral, BECAUSE it puts trapezoidal polygon under the curve. And for a convex function like sin(x), on an integrl where in theory the integral wants to be 2 exactly, trapz always misses by just a tiny bit, UNDER the true value.

So nothing you can do will allow trapz to solve the problem. You can never get all of the way there.

Now, suppose the goal was to find a point using numerical methids, such that the integral is 1, as close as trapz can do? That is doable.