marg_MAP option to use LBFGS hessians

[kimmywu]

This PR adds the option for marg_MAP to use LBFGS Hessians updates at each step, as well as terminating the MAP estimate either at nstep or when ϕtol is reached.

[kimmywu]

@marius311: For some unknown reason, converting to Map using Map(lastHg) is needed for the output lastHg to be scaled correctly when hess_method="lbfgs-hessian." Otherwise, it is orders of magnitude off (and looks like it's coming from a scale factor) when passed to history. It has the correct amplitude when applied to ϕ in the code.

[marius311]

Hmm I don't really remember tbh but glancing at

CMBLensing.jl/src/maximization.jl

Lines 200 to 212 in ec4bd10

	ϕ, = @⌛ optimize(
	objective,
	Map(ϕ),
	OptimKit.LBFGS(
	lbfgs_rank;
	maxiter = nsteps,
	verbosity = verbosity[1],
	linesearch = OptimKit.HagerZhangLineSearch(verbosity=verbosity[2], maxiter=5)
	);
	finalize!,
	inner = (_,ξ1,ξ2)->sum(unbatch(dot(ξ1,ξ2))),
	precondition = (_,η)->Map(Hϕ⁻¹*η),
	)

looks like I also have some Maps. I think in theory you could get rid of that by defining more of the things OptimKit needs like retract, inner (that one you already are), scale!, add!, and transport! as mentioned in their readme, although my guess is performance-wise the extra Map don't really matter so its probably fine.

[kimmywu]

I did see that you pass in Map-converted field variables in MAP_joint for the optimization in places. I tested both and in my case, they yield the same results with or without passing a Map-converted field and figure to just go without to reduce the back-and-forth.

Agreed that it doesn't slow down the code. But it does mean it's storing a larger vector (Map vs Fourier), and it's not the same (Fourier)type as the rest of the return keys. So I want to see if you already know of similar peculiar behavior. Or if this is a corner case, that I run into.