What are some ways to use trapz,integral or quad in reverse?
4 visualizzazioni (ultimi 30 giorni)
Mostra commenti meno recenti
Purush otham
il 30 Apr 2018
Modificato: Purush otham
il 2 Mag 2018
I have used ode45 to solve two simultaneous differential equations and plotted their solutions, such that : Area=trapz(x) or Area=integral(x), where 'x' is a function of 'phi'
My intention is to define the 'Area' first hand and find the associated 'phi'.
Example:
Let, Area of plot=trapz(x).
x=func(phi)
I give the Area and I want to find the associated 'phi'.
Another way to put it is:
If, V=trapz(x) for V=1, stop/print trapz(x)
0 Commenti
Risposta accettata
John D'Errico
il 30 Apr 2018
Modificato: John D'Errico
il 30 Apr 2018
It will essentially never happen that area is EXACTLY the value you want to find.
You could of course trivially use cumtrapz, to compute a cumulative trapezoidal integral. Then use find to find the first x where the area exceeds the target value.
x = linspace(0,2*pi,10);
y = cos(x);
Of course, we know that the integral is sin(x). But can we find the point where the integral is exactly 0.85?
Iy = cumtrapz(x,y)
Iy =
0 0.61647 0.94448 0.83056 0.32801 -0.32801 -0.83056 -0.94448 -0.61647 -4.4409e-16
It seems to lie between the second and third points, although there are several locations where it would happen. Find could find at least the first location where that value is exceeded.
ind = find(Iy > 0.85)
ind =
3
x(ind)
ans =
1.3963
Pretty crappy, really, in terms of the true value, which should arise at:
asin(0.85)
ans =
1.016
You point out that you understand tools like ODE45. Note that trapz and cumtrapz are equivalent to tools like Euler's method, so relatively poor ways to solve an ODE.
Could you write a cumulative trapezoidal rule that solves for the exact point where the trapezoidal integral hits some target value? Well, yes. It would still be terribly poor in terms of accuracy. Were I to want to do that (god knows why) I would probably write it as a piecewise linear spline. Then I would integrate the spline, and finally solve for the location(s) of interest.
But, if I was going to do it using a spline, just use a real cubic spline in the first place.
spl = spline(x,y);
Ispl = fnint(spl);
fnint comes with the curve fitting toolbox as I recall. If you lack that TB, it is easily enough written though.
Next, solve for the location of interest. In my case, I would just be lazy, and use the tools I've already written. Here, slmsolve, part of my SLM toolbox.
slmsolve(Ispl,0.85)
ans =
1.0126 2.1286
It found both locations where that event takes place.
But it is quite easy to solve using just fzero too.
target = 0.85;
fzero(@(x) ppval(Ispl,x) - target,1)
ans =
1.0126
fzero(@(x) ppval(Ispl,x) - target,2)
ans =
2.1286
Finally, could I have used an ODE solver, like ODE45 to solve this? Well, yes, but WHY would I? But yes, you could formulate an ODE to solve this problem using events.
Più risposte (0)
Vedere anche
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!