Depending on which one of n or p you want to eliminate in favor of sigma, you could use solve in combination with subs.

Here I'm eliminating n:

solve(Variance(X3)=sigma^2, {n}):
j3:=subs(%,k3);
 

Preben Alsholm

This method seems to work well.

I have minimized the residuals of the first 10 equations and treated the last two as equality constraints.

This gives an answer, which makes it possible for fsolve to succeed.

res2:=LSSolve(Res[1..10],{eq11,eq12,x1<=1,x2<=1,x3<=1,x4<=1,x5<=1,x6<=1,x7<=1,x8<=1,x9<=1,x10<=1,T>=Ti},
assume=nonnegative,iterationlimit=5000,feasibilitytolerance=.001,optimalitytolerance=0.001);
evalf(eval([eq1, eq2,eq3,eq4,eq5,eq6,eq7,eq8,eq9,eq10,eq11,eq12], res2[2]));
fsolve({eq1, eq2,eq3,eq4,eq5,eq6,eq7,eq8,eq9,eq10,eq11,eq12}, convert(res2[2],set));

{T = 3410.570845, n = 9.711378242, x1 = 0.03121886702, x10 = 0.01450927764,

  x2 = 0.1420646247, x3 = 0.6380365355, x4 = 0.009296381139,

  x5 = 0.07175312883, x6 = 0.04147335723, x7 = 0.01776536923,

  x8 = 0.006835808416, x9 = 0.02704665037}
 

Preben Alsholm

Here is my code after having executed your assignments of equations etc.

(lhs means 'left hand side: see ?lhs)

I have multiplied the last equation by 10^(-2) to reduce its weight in the optimization. This was the least weigt I had success with.

Preben Alsholm

with(Optimization):
Res:=map(lhs-rhs,[eq1, eq2,eq3,eq4,eq5,eq6,eq7,eq8,eq9,eq10,eq11,eq12*10^(-2)]):
#res:=LSSolve(Res,{x1<=1,x2<=1,x3<=1,x4<=1,x5<=1,x6<=1,x7<=1,x8<=1,x9<=1,x10<=1,T>=Ti},assume=nonnegative,iterationlimit=5000);
res:=LSSolve(Res,{x1<=1,x2<=1,x3<=1,x4<=1,x5<=1,x6<=1,x7<=1,x8<=1,x9<=1,x10<=1,T>=Ti},assume=nonnegative,iterationlimit=5000,optimalitytolerance=0.0001);
[0.00517923021408238112, [T = 2782.09143970567812, n = 10.1585637591277820,

  x1 = 0.00872658980230692216, x10 = 0.0748602181855385214,

  x2 = 0.0766666505398920090, x3 = 0.579521666596851892,

  x4 = 0.0000115210794556045040, x5 = 0.0901704459988488805,

  x6 = 0.0852556575745576917, x7 = 0.00577381343641316520,

  x8 = 0.00427475756541988076, x9 = 0.0642723480725257840]]
evalf(eval([eq1, eq2,eq3,eq4,eq5,eq6,eq7,eq8,eq9,eq10,eq11,eq12], res[2]));
[                                                              
[1. = 1.004651844, 4. = 4.001364580, 3.333333334 = 3.329171039,

  12.53333334 = 12.53468790, 0.9895336689 = 1., 0.02234128225 = 0.03162447862,

  0.0001704541989 = 0.02341381146, 0.001029281057 = 0.06427234807,

  0.0002380276970 = 0.07486021819, 0.07729879942 = 0.07666665054,

                                                5                 5]
  0.01208325906 = 0.008726589802, 6.693006871 10  = 6.693006871 10 ]
 

The underscore is often used by Maple when it has to come up with a global name, here a summation index.

In finding the general solution to a differential equation Maple creates names like _C1, _C2, etc.

I think the underscore is meant as a signal that this is a global variable invented by Maple, and don't use that kind of name yourself.

That the names are global implies that you can do e.g.

subs(_t=k,%);
                      infinity                        
                       -----                          
                        \     /   (1 - rho)          k\
                         )    |  k          q (1 - q) |
                        /     |- ---------------------|
                       -----  \        -1 + rho       /
                       k = 0                          
Preben Alsholm

C:=x->a(2*x^2+(V/x))+b(8*x+(V/x^2));
diff(C(x),x);
#The derivative as a function:
D(C);

# But maybe what you really mean is:

C:=x->a*(2*x^2+V/x)+b*(8*x+V/x^2);

Same syntax, though.

Preben Alsholm
 

You may try LSSolve in the Optimization package.

with(Optimization):
Res:=map(lhs-rhs,[eq1, eq2,eq3,eq4,eq5,eq6,eq7,eq8,eq9,eq10,eq11,eq12]):
res:=LSSolve(Res,{x1<=1,x2<=1,x3<=1,x4<=1,x5<=1,x6<=1,x7<=1,x8<=1,x9<=1,x10<=1,T>=Ti},assume=nonnegative,iterationlimit=1000);
Warning, limiting number of major iterations has been reached
  [199.025377293247600, [T = 2913.65979012579737, n = 1.89800495418788406,

    x1 = 0.419279311450406810, x10 = 0.489249736215187026,

    x2 = 0.403106848317799415, x3 = 0.985749025693950176,

    x4 = 0.330160746572949249, x5 = 0.242974425918598997,

    x6 = 0.189383244239871934, x7 = 0.372912480520804446,

    x8 = 0.353892475782923488, x9 = 0.393349253160286405]]

evalf(eval([eq1, eq2,eq3,eq4,eq5,eq6,eq7,eq8,eq9,eq10,eq11,eq12], res[2])):

fsolve({eq1, eq2,eq3,eq4,eq5,eq6,eq7,eq8,eq9,eq10,eq11,eq12}, convert(res[2],set));

The results of the last two lines show that you still don't have a solution to your equations, but at least the right and left hand sides are of the same order of magnitude.

Preben Alsholm

Here is a way to do it.

ode:=diff(y(x),x)=y(x);
sol1:=dsolve({ode,y(0)=1},numeric);
sol1(2);
#In this case op(2,sol1(x)):
Y1 := x-> rhs(op(2, sol1(x)));
Y1(.1234);
evalf(Int(Y1, 0 .. 1));
 

# or instead

sol2:=dsolve({ode,y(0)=1},numeric,output=listprocedure);
Y2:=subs(sol2,y(x));
evalf(Int(Y2, 0 .. 1));
 

My preferred method is the latter.

Preben Alsholm

The following also works, but is a lot slower than the method given by Robert Israel.

evalf(D[2](uu)(2,1));
plot(D[2](uu)(x,2),x=0..5,numpoints=5,adaptive =false);
evalf(Int(D[2](uu)(x,2),x=0..5,method=_Dexp,epsilon=1e-3));

Preben Alsholm
 

Generally speaking, I think this message is really bad news.

I guess you don't have a backup.

It is a good idea to save your document without output.  Remove output before saving.

Hopefully, someone can help you. But my experience tells me that there is not much to do about it.

Preben Alsholm

There is a very good reason why Maple cannot solve for the leading derivatives although they appear linearly in the system RedSys.

Consider the linear system with the unknowns being the leading derivatives, X16_1,X21_1,X63_2,X74_3,X97_2.

Then the coefficient matrix has only got rank 3. Had it been 5, there wouldn't have been a problem.

SYS:=subs(diff(X74(t),t,t,t)=X74_3,diff(X74(t),t,t)=X74_2,diff(X74(t),t)=X74_1,X74(t)=X74,diff(X97(t),t,t)=X97_2,diff(X97(t),t)=X97_1,X97(t)=X97,diff(X63(t),t,t)=X63_2,diff(X63(t),t)=X63_1,X63(t)=X63,diff(X21(t),t)=X21_1,X21(t)=X21,diff(X16(t),t)=X16_1,X16(t)=X16,RedSys):
L:=[op(1..4,SYS),lhs(SYS[5])=rhs5]:
A,b:=LinearAlgebra:-GenerateMatrix(L,[X16_1,X21_1,X63_2,X74_3,X97_2]):
LinearAlgebra:-Rank(A);
LinearAlgebra:-GaussianElimination(A);
LinearAlgebra:-Rank(<A|b>);

Preben Alsholm

Besides maybe other problems there is a problem with the initial conditions:

convert(RedSys,D):
eval(%,t=0):
res:=eval(%,incons1);
   [                                                                 
   [0 = 0.002, -800.0100000 @@(D, 2)(X97)(0) = 0, 0 = 0, 4.106328774 10  

                                                    
     @@(D, 2)(X63)(0) - 6.056424308 10   @@(D, 3)(X74)(0)

                                                      ]
      - 0.000002204071969 @@(D, 2)(X74)(0) = 0, 0. = 0]

Preben Alsholm

Another comment.

If you do

solve(RedSys,{diff(X74(t),t,t,t),diff(X97(t),t,t),diff(X63(t),t,t),diff(X21(t),t),diff(X16(t),t)});
 

NULL is returned.

Thus Maple cannot solve for the leading derivatives, here meaning the highest derivatives of the 5 functions appearing in the 5 equations.

Preben Alsholm

If you want to run muller.mw (not mulleralg.mw) as it is written, you need the library 'nalib',  which is available at the top of the page
http://www.cbu.edu/~wschrein/pages/M329Notes.html

Choose 'self-extracting Windows archive'. Save it in a convenient location. You need to know the path. Say it is

"F:/mymath/nalib"

then in Maple you do

restart;
libname:="F:/mymath/nalib", libname;
with(numanal);
#then you are ready.
Make sure that in the path you use `/` or  `\\`, not `\`.
 

The worksheet mulleralg.mw doesn't need the library because the procedure definition is right in there.

Walter Schreiner's numanal package contains many other procedures too.

Preben Alsholm

You may experiment with something like this.

f:=x^2: s:=1/2:
plot(s*f,x=0..4,scaling=constrained,tickmarks=[default,[seq(i=i/s,i=0..8)]]);
f:=x^2+y^3: s:=1/4:
plot3d(s*f,x=0..4,y=0..2,scaling=constrained,tickmarks=[default,default,[seq(i=i/s,i=0..6)]],axes=boxed);
 

Preben Alsholm

f:=Vector(5):
for i to 5 do
f[i] := 5.*i
end do:
f;
 

A:=Matrix(9,9,9);
sqrt~(A); #in Maple 13
map(sqrt,A); #in earlier versions, but also in 13.

B:=Matrix(9,symbol=b);
A/~B; #in Maple 13
zip(`/`,A,B); #in earlier versions, but also in 13.
 

Preben Alsholm