Clear["Global`*"]
SeedRandom[7];
n = 3; m =3;zMax = 1;
grid = Table[{x/n,y/m,zMax*Random[]},{x,0,n},{y,0,m}];
polygonalSurfacePlot=ListPlot3D[Flatten[grid,1],PlotStyle->{Opacity[0.25]},Mesh->None];
B[i_,u_,k_] := Binomial[k,i] u^i ( 1 - u)^(k-i);
p[u_,v_] = Sum[B[i,u,n]B[j,v,m]grid[[i+1,j+1]],{i,0,n},{j,0,m}];
grid[[n+1,m+1]];
p[1,1]
smoothPatchPlot=Plot3D[p[u,v][[3]],{u,0,1},{v,0,1},PlotStyle->{Blue,Opacity[0.75]},Mesh->None];
controlPointsPlot=ListPointPlot3D[Flatten[grid,1],PlotStyle->Directive[Orange,PointSize[0.03]]];
combined1=Show[polygonalSurfacePlot,smoothPatchPlot,controlPointsPlot];
e = Dot[D[p[u,v],u],D[p[u,v],u]];
g = Dot[D[p[u,v],v],D[p[u,v],v]];
f = Dot[D[p[u,v],u],D[p[u,v],v]];
ex1 = u2''[t] + D[e,u]/(2 e) + u2'[t]^2 + D[e,v]/e u2'[t] v2'[t] - D[g,u]/(2 e) v2'[t]^2/.{u->u[t], v->v[t],u2->u,v2->v}//Simplify;
ex2 = v2''[t] - D[e,v]/(2g) u2'[t]^2 + D[g,u]/g u2'[t] v2'[t] + D[g,v]/(2g) v2'[t]^2/.{u->u[t], v->v[t],u2->u,v2->v}//Simplify;

pInitial ={0.5,0.5}; (* These are u, v coordinates *)
pFinal = {0.75,0.15};
pInitialPatch = p[pInitial[[1]],pInitial[[2]]];
pFinalPatch = p[pFinal[[1]],pFinal[[2]]];
targetPtsPlot=ListPointPlot3D[{{pInitialPatch},{pFinalPatch}},PlotStyle->{Directive[Green,PointSize[0.03]],Directive[Red,PointSize[0.03]]},PlotRange->{{0,1},{0,1},{0,zMax}}];

Sol[vx0_,vy0_] := NDSolve[{ex1==0,ex2==0,u[0]==pInitial[[1]],v[0]==pInitial[[2]],u'[0]==vx0,v'[0]==vy0},{u,v},{t,0,1}]

uend[vx0_,vy0_] := (u[1]/.Sol[vx0,vy0][[1]])/;vx0 ∈ Reals && vy0 ∈ Reals;
vend[vx0_,vy0_] := (v[1]/.Sol[vx0,vy0][[1]])/;vx0 ∈ Reals && vy0 ∈ Reals;

vx0Guess =((pFinalPatch-pInitialPatch)/Norm[pFinalPatch-pInitialPatch])[[1]]
vy0Guess =((pFinalPatch-pInitialPatch)/Norm[pFinalPatch-pInitialPatch])[[2]]
roots=FindRoot[{uend[vx0,vy0]==pFinal[[1]], vend[vx0,vy0]==pFinal[[2]]},{{vx0,vx0Guess},{vy0,vy0Guess}},MaxIterations->250];
sol=Sol[vx0/.roots[[1]],vy0/.roots[[2]]];
{u[t_],v[t_]} = {u[t],v[t]}/.sol[[1]];
solutionCurvePlot = ParametricPlot3D[{u[t],v[t],p[u[t],v[t]][[3]]},{t,0,1},PlotStyle->Red];
combined2=Show[smoothPatchPlot,solutionCurvePlot,targetPtsPlot,ViewPoint->{1.2`,2.4`,1.5`}]

Export["/Users/canderson/bezier_soln2.png",%257,"PNG",ImageResolution->250]

Export["bezier_setup1.png",combined1,ImageResolution->250]

Export["bezier_soln2.png",combined2,ImageResolution->250]