Look at this question and the answer to it.

This is a well-known bug in some versions of Maple <= 2018. As a workaround, besides the depthlimit option, you need to increase  Digits :

restart;
kernelopts(version);
Digits:=20:
Optimization:-LPSolve(4*x1+5*x2+5*x3+2*x4,
                      {33*x1 + 49*x2 + 51*x3 + 22*x4 <= 120},
                      maximize, assume ={nonnegint}, depthlimit=10);

        Maple 2018.2, X86 64 WINDOWS, Oct 23 2018, Build ID 1356656             
                            [13, [x1 = 2, x2 = 1, x3 = 0, x4 = 0]]

It works:

Mat:=module()
    description "My Package";
    option package;
    uses Statistics;
    export LinModel; 
    LinModel :=proc(xliste::list, yliste::list,x::algebraic, $)
        return Fit(a*x+b,xliste,yliste,x)
        end proc; # LinModel 

    end module; #Mat

xliste:=[1,2,3]: yliste:=[2,4,5]:
Mat:-LinModel(xliste,yliste, x);
 

Only numerically:

f:=x->x^5 + 2*x^3:
g:=y->fsolve(y = x^5 + 2*x^3, x);  # This is the inverse function for f

Examples of use:

g(1);
plot([f, g], -3..3, -3..3, color=[red,blue], scaling=constrained);
evalf(Int(g, 0..2));
D(g)(2);
                                    

`/`(a, b); 
                       a/b

A simple procedure  LongSum  solves your problem. Formal parameters:  f  -  procedure for the common term, n0 - number of terms before the ellipsis,  n - position of the last term.

LongSum:=proc(f,n0,n)
`%+`(seq(f(k),k=1..n0))+`...`+f(n);
end proc:

Examples of use:
LongSum(n->x[n], 2, n);
LongSum(n->n/(n+1), 3, k);
LongSum(n->n/(n+1), 3, 100);

                           

If there is no any formula for the common term, then you can simply write (without a procedure) as in the example below:

a+b+c+`...`+f;

It is easy to modify this procedure to display not one, but several last terms or some intermediate terms, etc.

Download Check.mw

f2:=cos(y-Pi/3):
is(expand(f1/f2)=1);

                         true

Formally, the expressions  f1  and  f2  do not coincide for all values x and y. For example we get an undefined expression  0/0  for x=0, y=0 :
 

eval([numer(f1),denom(f1)],[x=0,y=0]);

                         [0, 0]
 

m := n->piecewise(type(n,even),0, 1):
f := proc(a,b)
if type(a*b,integer) then return m(a*b) else FAIL fi;
end:
 
f(a,b);
f(2,3);
f(1,3);
f(a,3);

sin(30*Pi/180);
30^`°`;

Example:

restart;
CurveFitting:-LeastSquares([[0, 1], [1, 2], [2, 3], [3, 10]], x);
A, B:=coeffs(%);

                               1/5+14*x*(1/5)
                             A, B := -1/5, 14/5

Good question - vote up.

restart;
P:=proc(n::nonnegint) 
local x, x0, v0, q, v:= <<x>>, s:= <<0>>, p:= <<1>>; 
   for x0 to n do v0:=eval(v,x=x0); p:= eval(p.v0,x=x0); s:= eval(s+v0, x=x0) od; 
   p, s; 
end proc:

P(6);

I have never used the command  rtable_eval  and do not know the reasons for this error.  A workaround below:

restart;
P:=proc()
local v;
v:= <x>;  
eval(v, x=1);
end proc():
P();

`&Delta;a`;
                                    

Should be  u*v  instead of  uv .

diff(BesselJ(alpha, sqrt(u^2+v^2-2*u*v*cos(phi))), u);