How to speed up code with the new compiler?

Above is the painting by Parmigiano “Self-portrait in a Convex Mirror”. Based upon my past contribution on reflection and anamorphism in spherical mirrors, I investigated if I could make a similar "painting"using Mathematica. Here is the geometrical setup:

C is the center of he mirror sphere, V is the view point and let S be a point to be reflected. The first problem is to find a point Q on the sphere such that the light ray VQ and the reflected ray QS form equal angles with the normal n on the sphere at Q.
I wrote a function including an NSolvevalues stating the laws of reflection: V,C, S and Q in same plane and vector angles VQ to n and SQ to n must be equal:

findIntersection[ptS : {xs_, ys_, zs_}, xv_] := 
 Module[{ptC = {0, 0, 0}, plane, ptV = {xv, 0, 0}}, 
  plane = InfinitePlane[{ptC, ptV, ptS}]; 
  Quiet@First[
    NSolveValues[{x, y, z} \[Element] Sphere[] && 
      VectorAngle[ptV - {x, y, z}, {x, y, z}] == 
       VectorAngle[ptS - {x, y, z}, {x, y, z}] && {x, y, 
        z} \[Element] plane && x > 0, {x, y, z}]]]

Works all right but the problem is the time it takes to find just one point Q. And I need a Q for every pixel in an image!

RepeatedTiming[findIntersection[{3.5, -4.5, -1.5}, 7.55]]

ScreenShot2022-01-03at11.20.06.png

I know about the new compiler of the Wolfram Language but I am not familiar enough with it to write myself a CompiledFunction version of my function. The old Compile does not work properly for this and takes even longer:

findIntersectionCF = Compile[{{ptS, _Real, 1}, {xv, _Real}},
  First@NSolveValues[{x, y, z} \[Element] Sphere[] && 
     VectorAngle[ptV - {x, y, z}, {x, y, z}] == 
      VectorAngle[ptS - {x, y, z}, {x, y, z}] && {x, y, z} \[Element] 
      InfinitePlane[{{0., 0., 0.}, {xv, 0., 0.}, ptS}] && x > 0, {x, 
     y, z}], CompilationTarget -> "C"]

5468ScreenShot2022-01-03at11.22.10.png

My question is : Can somebody at Wolfram or this community convert my function to the new CompiledFunction code and speed my function up?

I made an error in the formulation of my question. It should read: “Can somebody at Wolfram or this community convert my function to the new CompiledCodeFunction (instead of CompileFunction code)and speed my function up?” Sorry about that…

Not sure if the “new” compiler (FunctionCompile) can be of use here. The workhorse function (NSolveValues) is already a C level function in the Wolfram Language and thus as fast as can be.

It is much faster to avoid region specifications and inequalities when feasible to do so. Not sure how much it mattered below but I just made the sphere and plane conditions into explicit algebraic equations.

Also removing transcendentals by taking inverse functions, removing radicals by squaring, clearing denominators all can improve speed. At the end you impose the removed conditions e.g. with Cases and/or Select. Actually I found it expedient to use both.

findIntersection[ptS : {xs_, ys_, zs_}, xv_] := Module[
  {ptC = {0, 0, 0}, ptV = {xv, 0, 0}, angV, angS, solns, x, y, z},
  angV = Simplify[Cos[VectorAngle[ptV - {x, y, z}, {x, y, z}]]^2, 
    Element[{x, y, z}, Reals]];
  angS = Simplify[Cos[VectorAngle[ptS - {x, y, z}, {x, y, z}]]^2, 
    Element[{x, y, z}, Reals]];
  solns = 
   Cases[Quiet@
     NSolveValues[{x^2 + y^2 + z^2 == 1, 
       Numerator[Together[angV - angS]] == 0, zs*y - ys*z == 0}, {x, 
       y, z}], {x_Real, _Real, _Real} /; x > 0];
  SelectFirst[solns, (Abs[
       VectorAngle[ptV - #, #] - VectorAngle[ptS - #, #]] < 10^(-6)) &]
  ]

In[38]:= Timing[findIntersection[{3.5, -4.5, -1.5}, 7.55]]

(* Out[38]= {0.038289, {0.888266, -0.435758, -0.145253}} *)

Couple orders of magnitude.

Thank you Daniel, excellent and useful suggestions. There is something for me to learn in every step of your solution!

Thanks for your reply Arnaud. As suggested by Daniel, it appears better in some cases to rephrase the function than to Compile it.