Wednesday, March 24, 2021

Must be solved iteratively

One annoyance with h>3 is that Russell and Ocampo dumped a must be solved iteratively function on me. Find λ where f(λ) = 0. Ick!

But there was good news. f(λ) was differentiable! That meant I could call ol' Isaac to my rescue. I just had to differentiate that f(λ) monster on λ. And to turn both f(λ) and my new f'(λ) into anonymous functions of course. Anyway, for your perusal:

private delegate double SolveForL(double l);
private static double IterateNewtonRaphson(
SolveForL fun, 
SolveForL dfun, 
double startguess, 
double epsilon) {
 double l;
 double lNext = startguess;
 int iter = 0;
 do {
  l = lNext;
  lNext = l - fun(l) / dfun(l);
  iter++;
  if (iter > 1000) throw new ApplicationException("won't halt");
 } while (Math.Abs(lNext - l) > epsilon);
 return l;
}

My initial guess could fail, famously with sinusoidals; Mars/Earth 2-4-1-4 failed with initial zero. To that: more good news! These guys assured me that lambda would be bound, here between zero and right-angles. So I wrote this wrapper, to find (hopefully) the first real solution, and leave with it:

private static double IterateBoundedNewtonRaphson(
SolveForL fun, 
SolveForL dfun, 
double lower,
double upper,
double epsilon) {
 double l = 0; //dummy for now
 double guess = lower + double.Epsilon;
 do {
   try {
   l = IterateNewtonRaphson(fun, dfun, guess, epsilon);
   } catch(ApplicationException) {
                    //do nothing, it's a nonsolution
   }
   guess += 0.1;
 } while ((l <= lower || l > upper) && guess < upper);
 if (guess >= upper) throw new ApplicationException("No Solution");
 return l;
}

Just var inr = IterateBoundedNewtonRaphson(fun, dfun, 0, Math.PI / 2.0, epsilon); and we're off. 2-4-1-4 has lambda 0.9960 radians, turns out.

I did a sanity-check thingy after sending these functions over there, something like if (Math.Abs(eq8(inr)) > 10F * epsilon) throw new ApplicationException("you goofed");.

Still not saying I've totally sussed the h > 3 problem set. Certainly not saying the above is the efficient way to go about all of this. I wasn't even sure I was doing it right which is why I had that sanity-check.

But it is looking better.

No comments:

Post a Comment