You are trying to solve the ellipsoidal triangle where
φ1, α2, s12
are given. (See Figure 1 of the Wikipedia article, Geodesics on an
ellipsoid, for the notation.) In your case, we have
α2 = ±½π because of the requirement of
tangency. This is Problem 7 in §10 of my paper
Geodesics on an ellipsoid of revolution (Feb. 2011)
The solution entails assuming a value of α1 (perhaps by
solving an equivalent spherical problem), solving the ellipsoidal
triangle given
φ1, α1, α2
(Problem 2), determining the resulting s12, and correcting
α1 using Newton's method. The derivative needed for
Newton's method, ds12/dα1 is given by
Eq. (79). Problem 2 is trivial to solve: convert to the auxiliary
sphere (β1, α1, α2),
determine σ12 using some subset of Eqs. (9) thru (18),
and find s12 by using the routine Geodesic::ArcDirect in my
library GeographicLib. (Here β is the parametric latitude and
σ is the arc distance on the auxiliary sphere.)
Assuming that your initial guess is sufficiently good, this will
converge to the desired solution quadratically (the number of correct
digits in α1 will double on each iteration).