@Kitonum Error still there  

Error, (in simpl/reloprod) invalid terms in product: (0 = 121.4778145)^2

@Preben Alsholm 

Hi, I change to fsolve but it shows

fsolve({g1 = 0, g2 = 0, g3 = 0});
Error, (in fsolve) number of equations, 3, does not match number of variables, 4.

But I only have three variables there.

Thanks for your help!

@Preben Alsholm

Hi there,

I guess gam1 will be positive, but gam0 not for sure.

I put some initial values for n1,n2 and n3, 

so I change the code as:

s := 1/(273.16+50); s1 := 1/(273.16+145); s3 := 1/(273.16+250); s2 := 1/(273.16+197.5); gamma0 := 0.1e-3; gamma1 := .5; gamma2 := 0.15e-2; beta := -3800;
c := 300; n := 200; tau1 := 99; tau2 := 120;
Delta := solve(1-exp(-(gam0*tau1+(1/2)*gam1*tau1^2)*exp(beta*s1)) = 1-exp(-(gam0*a+(1/2)*gam1*a^2)*exp(beta*s2)), a);
a := Delta[1];
Theta := solve(1-exp(-(gam0*(a+tau2-tau1)+(1/2)*gam1*(a+tau2-tau1)^2)*exp(beta*s2)) = 1-exp(-(gam0*b+(1/2)*gam1*b^2)*exp(beta*s3)), b);
b := Theta[1];

lambda0 := gamma2*t^2+gamma1*t+gamma0;
n1 := evalf(n*(int(lambda0*exp(beta*s1)*exp(-(gamma0*t+(1/2)*gamma1*t^2+(1/3)*gamma2*t^3)*exp(beta*s1)), t = 0 .. tau1)));
56.49262074
n2 := evalf((n-n1)*(Int(lambda0*exp(beta*s2)*exp(-(gamma0*t+(1/2)*gamma1*t^2+(1/3)*gamma2*t^3)*exp(beta*s2)), t = a .. a+tau2-tau1)));
n3 := (n-n1-n2)*(Int(lambda0*exp(beta*s3)*exp(-(gamma0*t+(1/2)*gamma1*t^2+(1/3)*gamma2*t^3)*exp(beta*s3)), t = b .. c));
n4 := n-n1-n2-n3;

f1 := (gamma2*t^2+gamma1*t+gamma0)*exp(beta*s1)*exp(-(gamma0*t+(1/2)*gamma1*t^2+(1/3)*gamma2*t^3)*exp(beta*s1));
f2 := (gamma0+gamma1*(a+t-tau1)+gamma2*(a+t-tau1)^2)*exp(beta*s2)*exp(-(gamma0*(a+t-tau1)+(1/2)*gamma1*(a+t-tau1)^2+(1/3)*gamma2*(a+t-tau1)^3)*exp(beta*s2));
f3 := (gamma0+gamma1*(b+t-tau2)+gamma2*(b+t-tau2)^2)*exp(beta*s3)*exp(-(gamma0*(b+t-tau2)+(1/2)*gamma1*(b+t-tau2)^2+(1/3)*gamma2*(b+t-tau2)^3)*exp(beta*s3));

g1 := -n1*(Int((1/(gam1*t+gam0)-t*exp(theta*s1))*f1, t = 0 .. tau1))-n2*(Int((1/(gam0+gam1*(a+t-tau1))-(a+t-tau1)*exp(theta*s2))*f2, t = tau1 .. tau2))-n3*(Int((1/(gam0+gam1*(b+t-tau2))-(b+t-tau2)*exp(theta*s3))*f3, t = tau2 .. c))+(n-n1-n2-n3)*(1/(gam0+gam1*(b+c-tau2))-(b+c-tau2)*exp(theta*s3))*(gamma0+gamma1*(b+c-tau2)+gamma2*(b+c-tau2)^2)*exp(beta*s3)*exp(-(gamma0*(b+c-tau2)+(1/2)*gamma1*(b+c-tau2)^2+(1/3)*gamma2*(b+c-tau2)^3)*exp(beta*s3));

g2 := -n1*(Int((t/(gam1*t+gam0)-(1/2)*t^2*exp(theta*s1))*f1, t = 0 .. tau1))-n2*(Int(((a+t-tau1)/(gam0+gam1*(a+t-tau1))-(1/2)*(a+t-tau1)^2*exp(theta*s2))*f2, t = tau1 .. tau2))-n3*(Int(((b+t-tau2)/(gam0+gam1*(b+t-tau2))-(1/2)*(b+t-tau2)^2*exp(theta*s3))*f3, t = tau2 .. c))+(n-n1-n2-n3)*((b+c-tau2)/(gam0+gam1*(b+c-tau2))-(1/2)*(b+c-tau2)^2*exp(theta*s3))*(gamma0+gamma1*(b+c-tau2)+gamma2*(b+c-tau2)^2)*exp(beta*s3)*exp(-(gamma0*(b+c-tau2)+(1/2)*gamma1*(b+c-tau2)^2+(1/3)*gamma2*(b+c-tau2)^3)*exp(beta*s3));

g3 := -n1*s1+n1*(Int(s1*(gam0*t+(1/2)*gam1*t^2)*exp(theta*s1)*f1, t = 0 .. tau1))-n2*s2+n2*(Int(s2*(gam0*(a+t-tau1)+(1/2)*gam1*(a+t-tau1)^2)*exp(theta*s2)*f2, t = tau1 .. tau2))-n3*s3+n3*(Int(s3*(gam0*(b+t-tau2)+(1/2)*gam1*(b+t-tau2)^2)*exp(theta*s3)*f3, t = tau2 .. c))+(n-n1-n2-n3)*s3*(gam0*(b+c-tau2)+(1/2)*gam1*c^2)*exp(theta*s3)*(gamma0+gamma1*(b+c-tau2)+gamma2*(b+c-tau2)^2)*exp(beta*s3)*exp(-(gamma0*(b+c-tau2)+(1/2)*gamma1*(b+c-tau2)^2+(1/3)*gamma2*(b+c-tau2)^3)*exp(beta*s3));

 solve({g1 = 0, g2 = 0, g3 = 0}, {gam0, gam1, theta})

and it shows like this,

Error, (in solve) cannot solve expressions with Int((7006148439/10000000000000)*((3/2000)*t^2+(1/2)*t+1/10000)*exp(-(7006148439/100000000000000000)*t-(175153711/1000000000000)*t^2-(175153711/500000000000000)*t^3), t = (3/2500000000000000000)*(-833333333200000000*gam0+(713658873/1250000000000000000)*(2130473782000000000000000000000000000000000000000000000*gam0^2+3154676655000000000000000000000000000000*gam0*(970987338*gam0^2+69767150990*gam0*gam1+3453473973000*gam1^2)^(1/2)+127716247600000000000000000000000000000

did I do something wrong there?

@Preben Alsholm 

Hi,

why you know gam0=0.01,gam1=3.62

evalf(eval([g1,g2],{gam0=0.01,gam1=3.62}));

and  you had no success using fsolve as in
fsolve({g1=0,g2=0}, {gam0=0..230, gam1=0.1..10});

is it possible that the solution is not in this range?

@Carl Love 

Oh! I made a mistake there for n3

it should be n3 := eval(n33, b1 = b); rather than n3 := eval(n33, a1 = a);

a and b both are the function of gam0 and gam1, so I don't have numeric values for it...

Is there any method can help me to solve it?

Furthermore, you said it's a very bad coding practice to use t both as a free variable and a variable of integration. I can not really understand this, could you explain it by example, please. Thanks !

@Carl Love 

I use

minimize(abs(int(0.1e-3+.5*t+0.2e-2*t^2-b*t-a, t = 0 .. 300)), location = true);

and the result is:

0., {[{a = 135.0001000 - 150. b, b = b}, 0.]}

stiil cannot get the numerical answer of b?

@Carl Love 

how to add the constraints in maple code if i wanna minimize the absolute value of int(0.1e-3+.5*t+0.2e-2*t^2-b*t-a, t = 0 .. 300)

@Preben Alsholm 

Yes, I got the same answer as what you got.

@Carl Love 

yeah, I have changed 

NLPSolve(1/(n^3*(F0*F1*F2-F0*F12^2-F01^2*F2+2*F01*F02*F12-F02^2*F1)), {u1 >= 40, u2 >= 20});

But the error still comes:

Error, (in Optimization:-NLPSolve) invalid arguments

@acer 

Hi acer, I try this way. But sometimes the answer doesn't seem right. Cause I have three variables, I can not use what the 3D plot looks like over that range to consider. Is there another way can solve optimize?

What if I want to add constrant in NLPSolve?

I set  u1 := n*(int(lambda0(t)*exp(beta*s1)*exp(-(gamma0*t+(1/2)*gamma1*t^2+(1/3)*gamma2*t^3)*exp(beta*s1)), t = 0 .. tau1))

u2 := (n-u1)*(int(lambda0(t)*exp(beta*s2)*exp(-(gamma0*t+(1/2)*gamma1*t^2+(1/3)*gamma2*t^3)*exp(beta*s2)), t = a .. a+tau2-tau1))

NLPSolve(subs(int=Int,1/(n^3*(F0*F1-F1))),{u1>40,u2>30}, variables=[tau1,tau2,t2]);

Error, (in Optimization:-NLPSolve) strict inequalities are not valid constraints

@acer 

Hi acer,

I need to have tau1 tau2 and t2 as varibles to get there optimal values.

So I remove t2 in the beginning, and run the code again.

However, it doesn't work.

NLPSolve(1/(n^3*(F0*F1-F1)), tau1 = 115 .. 201, tau2 = 237 .. 273, t2=145..250);

How can I solve this problem?

@acer 

Hi, why we need to subs(int=Int,1/(n^3*(F0*F1-F1))

@Carl Love 

Thanks, it does work.

I also meet the same problem in the case below:

restart; with(linalg); with(stats); with(plots); with(Statistics); with(LinearAlgebra); with(Optimization);

lambda0 := proc (t) options operator, arrow; gamma0+gamma1*t+gamma2*t^2;
lambda := lambda0(t)*exp(beta*s);
s := 1/(273.16+50); s1 := 1/(273.16+145); s3 := 1/(273.16+250); s2 := 1/(273.16+(1/2)*(145+250)); gamma0 := 0.1e-3; gamma1 := .5; gamma2 := 0; beta := -3800;
c := 300; n := 200;

Theta := solve(1-exp(-(gamma0*tau1+(1/2)*gamma1*tau1^2)*exp(beta*s1)) = 1-exp(-(gamma0*a+(1/2)*gamma1*a^2)*exp(beta*s2)), a);
a := Theta[2];

Lambda := solve(1-exp(-(gamma0*(a+tau2-tau1)+(1/2)*gamma1*(a+tau2-tau1)^2)*exp(beta*s2)) = 1-exp(-(gamma0*b+(1/2)*gamma1*b^2)*exp(beta*s3)), b);

b := Lambda[1];

F0 := int(exp(beta*s1)*exp(-(gamma0*t+(1/2)*gamma1*t^2+(1/3)*gamma2*t^3)*exp(beta*s1))/(gamma2*t^2+gamma1*t+gamma0), t = 0 .. tau1)+evalf(int(exp(beta*s2)*exp(-(gamma0*(a+t-tau1)+(1/2)*gamma1*(a+t-tau1)^2+(1/3)*gamma2*(a+t-tau1)^3)*exp(beta*s2))/(gamma0+gamma1*(a+t-tau1)+gamma2*(a+t-tau1)^2), t = tau1 .. tau2))+evalf(int(exp(beta*s3)*exp(-(gamma0*(b+t-tau2)+(1/2)*gamma1*(b+t-tau2)^2+(1/3)*gamma2*(b+t-tau2)^3)*exp(beta*s3))/(gamma0+gamma1*(b+t-tau2)+gamma2*(b+t-tau2)^2), t = tau2 .. c));

F1 := int(t^2*exp(beta*s1)*exp(-(gamma0*t+(1/2)*gamma1*t^2+(1/3)*gamma2*t^3)*exp(beta*s1))/(gamma2*t^2+gamma1*t+gamma0), t = 0 .. tau1)+evalf(int((a+t-tau1)^2*exp(beta*s2)*exp(-(gamma0*(a+t-tau1)+(1/2)*gamma1*(a+t-tau1)^2+(1/3)*gamma2*(a+t-tau1)^3)*exp(beta*s2))/(gamma0+gamma1*(a+t-tau1)+gamma2*(a+t-tau1)^2), t = tau1 .. tau2))+evalf(int((b+t-tau2)^2*exp(beta*s3)*exp(-(gamma0*(b+t-tau2)+(1/2)*gamma1*(b+t-tau2)^2+(1/3)*gamma2*(b+t-tau2)^3)*exp(beta*s3))/(gamma0+gamma1*(b+t-tau2)+gamma2*(b+t-tau2)^2), t = tau2 .. c));
Warning, computation interrupted

F2 := int(t^4*exp(beta*s1)*exp(-(gamma0*t+(1/2)*gamma1*t^2+(1/3)*gamma2*t^3)*exp(beta*s1))/(gamma2*t^2+gamma1*t+gamma0), t = 0 .. tau1)+evalf(int((a+t-tau1)^4*exp(beta*s2)*exp(-(gamma0*(a+t-tau1)+(1/2)*gamma1*(a+t-tau1)^2+(1/3)*gamma2*(a+t-tau1)^3)*exp(beta*s2))/(gamma0+gamma1*(a+t-tau1)+gamma2*(a+t-tau1)^2), t = tau1 .. tau2))+evalf(int((b+t-tau2)^4*exp(beta*s3)*exp(-(gamma0*(b+t-tau2)+(1/2)*gamma1*(b+t-tau2)^2+(1/3)*gamma2*(b+t-tau2)^3)*exp(beta*s3))/(gamma0+gamma1*(b+t-tau2)+gamma2*(b+t-tau2)^2), t = tau2 .. c));
Warning, unable to determine if -1/5000 is between 0 and tau1; try to use assumptions or use the AllSolutions option

F01 := int(t*exp(beta*s1)*exp(-(gamma0*t+(1/2)*gamma1*t^2+(1/3)*gamma2*t^3)*exp(beta*s1))/(gamma2*t^2+gamma1*t+gamma0), t = 0 .. tau1)+evalf(int((a+t-tau1)*exp(beta*s2)*exp(-(gamma0*(a+t-tau1)+(1/2)*gamma1*(a+t-tau1)^2+(1/3)*gamma2*(a+t-tau1)^3)*exp(beta*s2))/(gamma0+gamma1*(a+t-tau1)+gamma2*(a+t-tau1)^2), t = tau1 .. tau2))+evalf(int((b+t-tau2)*exp(beta*s3)*exp(-(gamma0*(b+t-tau2)+(1/2)*gamma1*(b+t-tau2)^2+(1/3)*gamma2*(b+t-tau2)^3)*exp(beta*s3))/(gamma0+gamma1*(b+t-tau2)+gamma2*(b+t-tau2)^2), t = tau2 .. c));

NLPSolve(1/(n^3*(F0*F1*F2-F0+F01)), tau1 = 115 .. 201, tau2 = 237 .. 273);

In this case, i didn't use A,B,C,d,E function. 

I need to find tau1 and tau2 which can Minimize(1/(n^3*(F0*F1*F2-F0+F01)),,tau1=115..201, tau2 = 237..273})

But this error keeps coming :

Error, (in Optimization:-NLPSolve) integration range or variable must be specified in the second argument, got HFloat(1.0) = 121.0828419 .. HFloat(193.0828419)

@Carl Love 

restart; with(linalg); with(stats); with(plots); with(Statistics); with(LinearAlgebra); with(Optimization);
lambda0(t):=gamma0+gamma1*t:
lambda := lambda0(t)*exp(beta*s);
t1 := 145; t3 := 250; t2 := (t1+t3)*(1/2);
s := 1/(273.16+50); s1 := 1/(273.16+t1); s3 := 1/(273.16+t3); s2 := 1/(273.16+t2); gamma0 := 0.1e-3; gamma1 := .5; beta := -3800;

c := 300; n := 200; tau1 := 201;
Theta := solve(1-exp(-(gamma0*tau1+(1/2)*gamma1*tau1^2)*exp(beta*s1)) = 1-exp(-(gamma0*a+(1/2)*gamma1*a^2)*exp(beta*s2)), a);
a := Theta[2];

Lambda := solve(1-exp(-(gamma0*(a+tau2-tau1)+(1/2)*gamma1*(a+tau2-tau1)^2)*exp(beta*s2)) = 1-exp(-(gamma0*b+(1/2)*gamma1*b^2)*exp(beta*s3)), b);

b := Lambda[1];

A(c,n,m):=evalf(int(1/(y)*exp(-c*y^(2)),y=n..m))

B(c,n,m):=evalf(int(exp(-c*y^(2)),y=n..m))

C(c,n,m):=evalf(int(y*exp(-c*y^(2)),y=n..m))

d(c,n,m):=evalf(int(y^(2)*exp(-c*y^(2)),y=n..m))

E(c,n,m):=evalf(int(y^(3)*exp(-c*y^(2)),y=n..m))

F0 := exp(beta*s1)*exp(gamma0^2*exp(beta*s1)/(2*gamma1))*A((1/2)*gamma1*exp(beta*s1), gamma0/gamma1,gamma0/gamma1+tau1)/gamma1+exp(beta*s2)*exp(gamma0^2*exp(beta*s2)/(2*gamma1))*A((1/2)*gamma1*exp(beta*s2), gamma0/gamma1+a, gamma0/gamma1+tau2-tau1+a)/gamma1+exp(beta*s3)*exp(gamma0^2*exp(beta*s3)/(2*gamma1))*A((1/2)*gamma1*exp(beta*s3), gamma0/gamma1+b, gamma0/gamma1+c-tau2+b)/gamma1

F1 := exp(beta*s1)*exp(gamma0^2*exp(beta*s1)/(2*gamma1))*(gamma0^2*A((1/2)*gamma1*exp(beta*s1), gamma0/gamma1, gamma0/gamma1+tau1)/gamma1^2-2*gamma0*B((1/2)*gamma1*exp(beta*s1), gamma0/gamma1, gamma0/gamma1+tau1)/gamma1+C((1/2)*gamma1*exp(beta*s1), gamma0/gamma1, gamma0/gamma1+tau1))/gamma1+exp(beta*s2)*exp(gamma0^2*exp(beta*s2)/(2*gamma1))*(gamma0^2*A((1/2)*gamma1*exp(beta*s2), gamma0/gamma1+a, gamma0/gamma1+tau2-tau1+a)/gamma1^2-2*gamma0*B((1/2)*gamma1*exp(beta*s2), gamma0/gamma1+a, gamma0/gamma1+tau2-tau1+a)/gamma1+C((1/2)*gamma1*exp(beta*s2), gamma0/gamma1+a, gamma0/gamma1+tau2-tau1+a))/gamma1+exp(beta*s3)*exp(gamma0^2*exp(beta*s3)/(2*gamma1))*(gamma0^2*A((1/2)*gamma1*exp(beta*s3), gamma0/gamma1+b, gamma0/gamma1+c-tau2+b)/gamma1^2-2*gamma0*B((1/2)*gamma1*exp(beta*s3), gamma0/gamma1+b, gamma0/gamma1+c-tau2+b)/gamma1+C((1/2)*gamma1*exp(beta*s3), gamma0/gamma1+b, gamma0/gamma1+c-tau2+b))/gamma1

F01 := exp(beta*s1)*exp(gamma0^2*exp(beta*s1)/(2*gamma1))*(B((1/2)*gamma1*exp(beta*s1), gamma0/gamma1, gamma0/gamma1+tau1)-gamma0*A((1/2)*gamma1*exp(beta*s1), gamma0/gamma1, gamma0/gamma1+tau1)/gamma1)/gamma1+exp(beta*s2)*exp(gamma0^2*exp(beta*s2)/(2*gamma1))*(B((1/2)*gamma1*exp(beta*s2), gamma0/gamma1+a, gamma0/gamma1+tau2-tau1+a)-gamma0*A((1/2)*gamma1*exp(beta*s2), gamma0/gamma1+a, gamma0/gamma1+tau2-tau1+a)/gamma1)/gamma1+exp(beta*s3)*exp(gamma0^2*exp(beta*s3)/(2*gamma1))*(B((1/2)*gamma1*exp(beta*s3), gamma0/gamma1+b, gamma0/gamma1+c-tau2+b)-gamma0*A((1/2)*gamma1*exp(beta*s3), gamma0/gamma1+b, gamma0/gamma1+c-tau2+b)/gamma1)/gamma1

`Fβ` := int((s1^2*(gamma0*t+(1/2)*gamma1*t^2)*exp(beta*s1)*(gamma1*t+gamma0)*exp(beta*s1))*exp(-(gamma0*t+(1/2)*gamma1*t^2)*exp(beta*s1)), t = 0 .. tau1)+int((s2^2*(gamma0*(a+t-tau1)+(1/2)*gamma1*(a+t-tau1)^2)*exp(beta*s2)*(gamma0+gamma1*(a+t-tau1))*exp(beta*s2))*exp(-(gamma0*(a+t-tau1)+(1/2)*gamma1*(a+t-tau1)^2)*exp(beta*s2)), t = tau1 .. tau2)+int((s3^2*(gamma0*(b+t-tau2)+(1/2)*gamma1*(b+t-tau2)^2)*exp(beta*s3)*(gamma0+gamma1*(b+t-tau2))*exp(beta*s3))*exp(-(gamma0*(b+t-tau2)+(1/2)*gamma1*(b+t-tau2)^2)*exp(beta*s3)), t = tau2 .. c)+int((s3^2*(gamma0*(b+t-tau2)+(1/2)*gamma1*(b+t-tau2)^2)*exp(beta*s3)*(gamma0+gamma1*(b+t-tau2))*exp(beta*s3))*exp(-(gamma0*(b+t-tau2)+(1/2)*gamma1*(b+t-tau2)^2)*exp(beta*s3)), t = c .. infinity)

I need to have tau2 as varibles to get there optimal values ..

Minimize(1/((F0*F1-F01^2)*n^3*`Fβ`), tau2 = 237..273})