We need to find a point on the curve, in which  the normal passes through original point. The following code solves the problem for the original point  M=[0,-2] :

restart;

f:=x->1/x:

Eq:=(y=f(x0)-1/D(f)(x0)*(x-x0)):

M:=[0,-2]:

x0:=fsolve(eval(Eq, [x,y]=~M), x0=0..infinity);

A:=plot([1/x, f(x0)+D(f)(x0)*(x-x0), [t,eval(rhs(Eq),x=t), t=0..x0]], x=-2..4, -3..3, color=[red,blue,green], thickness=[2,1,2], discont):

B:=seq(plottools[disk](p,0.03,color=yellow), p=[M,[x0,f(x0)]]):

C:=plots[textplot]([[0.2,-2,"M"],[x0+0.2,f(x0)+0.1,P]], font=[TIMES,ROMAN,18]):

plots[display](A,B,C);  # Visualization - the point P is the projection of the point M onto the curve  y=1/x, x>0

 Addition: In general, the solution is not unique. But in this example, for each of the branches of the hyperbola (due to the convexity) solution is unique.

Acer method probably the shortest. Here are more 3 ways:

X := [1, 3, 8, 9, 10]:

`+`(seq(X[i]*ln(X[i]), i=1..nops(X)));

# Or

`+`(op(X*~(ln~(X)))); 

# Or in for loop

S := 0: for x in X do S := S+x*ln(x) end do: S;

                         21*ln(3)+24*ln(2)+10*ln(10)

                         21*ln(3)+24*ln(2)+10*ln(10)

                         21*ln(3)+24*ln(2)+10*ln(10)

To see the point of intersection, you can simply reduce the ranges of the axes at standard numpoints:

restart;

convert([92*x/790+18*x*(1-y)/1000=0.125,ln(46.59/x)=1/1.5*(ln(0.553/y)+2.5*ln((1-y)/(1-0.553)))], rational);

fsolve(%);

plots:-implicitplot(%%, x=0..1, y=0..0.05, thickness=2, color=[blue,red]);

Your system can easily be solved exactly. In general case, it will have  2^3 = 8  solutions. Since the orthogonal matrix  RT has determinant  1, the matrix  mat  should also be orthogonal with determinant  1. For solution we need to select 3 independent elements of the orthogonal matrix  mat (I took  mat[1,2], mat[1,3], mat[2,3]). Here are all the solutions for the specific orthogonal matrix  mat  with full check of one solution:

restart;

Rx := Matrix(3,3, [1,0,0,0,cos(a),-sin(a),0,sin(a),cos(a)]):   

Ry := Matrix(3,3, [cos(b), 0, -sin(b), 0, 1, 0, sin(b), 0, cos(b)]):  # Corrected matrix

Rz := Matrix(3,3, [cos(c), -sin(c), 0, sin(c), cos(c), 0, 0, 0, 1]):

RT := Rx.Rz.Ry:

mat:=<1/3,2/3,-2/3; -2/3,2/3,1/3; 2/3,1/3,2/3>;  # The specific orthogonal matrix with det=1

assign(seq(seq(g[i,j]=(RT[i,j]=mat[i,j]), j=1..3), i=1..3));

Sol:=allvalues([solve({g[1,2], g[1,3], g[2,3]})]);  # Found solutions

eval([seq(seq(g[i,j], j=1..3), i=1..3)], Sol[1,1]);  # The check of  Sol[1,1]  for all equations of the system

Additions.  1) You can verify that  Sol[4,1]  is also a solution. The rest of the solutions are not solutions of the whole system.

2) The code was edited.

Your procedure works for me, but the same can be written a much shorter.

Your way:

dList:=proc(a::list, b::list)

    local res::list, i::integer, n::integer;

    description "add list elements and diff";

    n:=min(nops(a), nops(b)) ;

res:=[seq(0, i=1..n)];

for i from 1 by 1 to n do

     res[i]:=diff(a[i]+b[i],x)

end do;

res;

end proc:

list1:= [x^3 , 2*x , 3]:

list2:= [x^5, 2*x^4 , 2*x]:

dList(list1,list2);

                    [3*x^2+5*x^4, 2+8*x^3, 2]

The shorter method (without user's procedure):

list1:= [x^3 , 2*x , 3]:

list2:= [x^5, 2*x^4 , 2*x]:

map(diff, list1+list2, x);

                     [3*x^2+5*x^4, 2+8*x^3, 2]

 Maybe your error is caused by previous assignments. Make restart before your procedure.

Everything to be beautiful, I changed all the colors and ranges for the axes. Instead of  Student package I just wrote the simple equations for the tangent and normal:

f := x-> x^2-2*x:

x0 := 1.5:

T := (D(f))(x0)*(x-x0)+f(x0):  # The tangant

N := -(x-x0)/(D(f))(x0)+f(x0):  # The normal

P := plottools[point]([x0, f(x0)], symbol = circle, symbolsize = 15):  # The point

Curves := plot([f(x), T, N], x = -4 .. 5, -4 .. 5, color = [red, blue, green], thickness = 2, legend = [curve, tangent, normal]):

plots[display](Curves, P, scaling = constrained);

You will have more opportunities to vary the shape of the knot, if you use it's parametric equations and  plots[tubeplot]  command. You are setting the parameters  rho>1  and coprime integers  p  and  q .

Example:

restart;

r := cos(q*phi)+rho:  x := r*cos(p*phi):  y := r*sin(p*phi):  z := -sin(q*phi):  p := 3:  q := 7:  rho := 3:

plots[tubeplot]([x, y, z], phi = 0 .. 2*Pi, radius = 0.3, style = surface, numpoints = 3000, scaling = constrained, axes = none);

See Carl's idea using  viewpoint  option  here

In this thread are also discussed some other ways to solve the problem.

There is substantial error in procedure  Simpson  - instead  n  should be  n/2 . I made a few less significant corrections. Now everything is working correctly:

Simpson := proc(f, a, b, n::even)

local h, S1, S2, S, i;

h := (b-a)/n;

S1 := 0.;

for i from 0 to n/2-1 do

S1 := S1 + f(a + (2*i+1)*h);

end do;

S2 := 0.;

for i from 1 to n/2-1 do

S2 := S2 + f(a + (2*i)*h);

end do;

evalf(h/3 * ( f(a)+f(b) + 4*S1 + 2*S2 )); 

end proc:

Example of work of the procedure  Simpson  and independent calculation:

Digits := 5;                          

f := x -> 1/sqrt(39.24*x-44.65*(x*arccos(x)-sqrt(1-x^2)-13.88*(1-x^2)^1.5));

Simpson(f, 0, 1, 100);

evalf(Int(f(x), x=0..1));

is(diff(abs(x), x) = x/abs(x))  assuming x < 0;

is(diff(abs(x), x) = x/abs(x))  assuming x > 0;

                                  true

                                  true

Slightly more programmatically:

op(-1, eval(splcurve) assuming v < 4.2);

If you want to get a normal manual solution, you can use  Student[Calculus1][ShowSolution]  command.

Example:

Student[Calculus1][ShowSolution](Int(sin(x)*(x + sin(x)),x));

roll:=rand(1..6):

L:=[seq([seq(roll(), j=1..10)], i=1..30)];

NumbersOfThrees:=map(t->ListTools[Occurrences](3, t), L);

The number  pi  should be coded as  Pi .

All curves are in the same plot:

restart;

Eq:=ln(2*s*t+2*s*x*sqrt(t)/sqrt(Pi)+x^2)-2*s*(arctan((s+x*sqrt(Pi/t))/sqrt(2*Pi*s-s^2))-(1/2)*Pi)/sqrt(2*Pi*s-s^2) = 0;

S:=[0.01,0.005,0.003,0.002];

plots[implicitplot]([seq(eval(Eq,s=S[i]), i=1..nops(S))], t=0..5, x=-2..2, color=[red,blue,green,yellow], gridrefine=5);

It is easy to verify identity

is(sin(2*x)^2 - cos(8*x)^2 = - cos(10*x)*cos(6*x));

                                         true

Therefore, the equation splits into two equations, which are easy to solve

solve({x >= 0, `or`(cos(10*x) = 0, cos(6*x)+1/2 = 0), x <= (1/2)*Pi}, AllSolutions, explicit);