vv

13827 Reputation

20 Badges

9 years, 316 days

MaplePrimes Activity


These are answers submitted by vv


 

restart;

 

The exact (symbolic) absolute min and max.
Two variables. The constraints must define a compact domain.
Based on extrema and solve (which are supposed to work).

absminmax:=proc(f::algebraic, constr::set(relation),candi::name:='_candidates')
local bdr, n,P,i,j,u,s,cand,x,y;
s:=indets([f,constr],name) minus {constants};
if nops(s)<>2 then error "Two indeterminates only" fi;
bdr:=seq(lhs(u)-rhs(u)=0,u=constr):
x,y:=s[];  n:=0:  P:=Vector():
for j to nops(constr) do for i to j-1 do
  s:=solve([bdr[i],bdr[j]],{x,y},explicit);
  for u in [s] do n:=n+1; P(n):=u od;
od od:
for i to nops(constr) do
  extrema(f, {bdr[i]},{x,y}, 's'):
  for u in s do n:=n+1: P(n):=u od;
od:
s:=solve([diff(f,x),diff(f,y)],{x,y},explicit):
for u in [s] do n:=n+1: P(n):=u od:
cand:=select(u -> type(rhs~(u),set(realcons)), {entries(P,nolist)});
cand:=select(u -> is( `and`( eval(constr, u)[]) ), cand);
candi:=cand; if cand={} then return {} fi;
min[index]( [seq(eval(f,u),u=cand)]):  
  s:=`min`=eval(f,cand[%]), cand[%];
max[index]( [seq(eval(f,u),u=cand)]):  
  s, `max`=eval(f,cand[%]), cand[%];
end:

 

Examples

 

 

f1:=2*x^2+4*x*y+2*x-2*y^2+y+4:

constr1:={ y>=(1/21)*x-283/42,
          x<=(3/29)*y+329/58,
          y<=(1/9)*x+131/18,
          x>=-(1/9)*y-113/18 }:

absminmax(f1, constr1);

min = -137657/784, {x = -2267/392, y = 2601/392}, max = 379/2, {x = 13/2, y = 8}

(1)

_candidates;

{{x = -7, y = 13/2}, {x = 5, y = -13/2}, {x = -2267/392, y = -1745/392}, {x = -2267/392, y = 2601/392}, {x = -11/2, y = -7}, {x = -3/8, y = -1/8}, {x = 13/2, y = 8}}

(2)

f2:=4*x*y-2*y^2+4*x+y+4:
constr2:={y >= -(1/13)*x-85/13, y <= -(25/14)*x+26/7, y <= (23/12)*x+89/12}:

absminmax(f2, constr2,'ccc');

min = -245, {x = 6, y = -7}, max = 66, {x = -7, y = -6}

(3)

absminmax(x^2+y^2, {x<=1,y<=3,x>=0,y>=0});

min = 0, {x = 0, y = 0}, max = 10, {x = 1, y = 3}

(4)

absminmax(x+y, {(x-1)^2+y^2<=1, y>=x^2});

min = 0, {x = 0, y = 0}, max = 2, {x = 1, y = 1}

(5)

absminmax(exp(x+y), {(x-1)^2+y^2<=1, y>=x^2});

min = 1, {x = 0, y = 0}, max = exp(2), {x = 1, y = 1}

(6)

absminmax(x+y, {(x-1)^2+y^2<=1, y>=x^2+10});

{}

(7)

absminmax(x^2-y^2, {(x-1)^2+y^2<=1, y>=x^2});

min = -1/2, {x = 1/2, y = (1/2)*3^(1/2)}, max = 1/4, {x = (1/2)*2^(1/2), y = 1/2}

(8)

absminmax(x^4-y^2, {(x-1)^2+y^2<=1, y<=x^2});

min = ((1/6)*(54+6*87^(1/2))^(1/3)-1/(54+6*87^(1/2))^(1/3))^4-(1/36)*(-(54+6*87^(1/2))^(4/3)+12*(54+6*87^(1/2))^(2/3)-72*(54+6*87^(1/2))^(1/3)+72*87^(1/2)+612)/(54+6*87^(1/2))^(2/3), {x = (1/6)*(54+6*87^(1/2))^(1/3)-1/(54+6*87^(1/2))^(1/3), y = -(1/6)*((54+6*87^(1/2))^(2/3)*(-(54+6*87^(1/2))^(4/3)+12*(54+6*87^(1/2))^(2/3)-72*(54+6*87^(1/2))^(1/3)+72*87^(1/2)+612))^(1/2)/(54+6*87^(1/2))^(2/3)}, max = 16, {x = 2, y = 0}

(9)

 


Download absminmax-sent.mw

Diff(int(value(z), t), t);


 

restart;

f:=((z^4*exp(2*z)+1)/(z+I)^3) - ((z^3+z)/((z-2*I)*(z-5))) + 8*Pi*exp(1);  # probably

(z^4*exp(2*z)+1)/(z+I)^3-(z^3+z)/((z-2*I)*(z-5))+8*Pi*exp(1)

(1)

poles:=[-I, 2*I, 5];

[-I, 2*I, 5]

(2)

C:= z -> is( abs(z-1) < sqrt(11/2) );

proc (z) options operator, arrow; is(abs(z-1) < sqrt(11/2)) end proc

(3)

2*Pi*I*add( residue(f, z=a), a=select(C, poles) );

(2*I)*Pi*(12/29-(30/29)*I+(102+112*I+(-2-4*I)*((4*I)*(exp(I))^2+5*(exp(I))^2)/(exp(I))^2)/(-(15*I)*(exp(I))^2+3*(exp(I))^2))

(4)

evalf(%);

4.564531460+58.76512338*I

(5)

 


 

The fact is that pdsolve is not very competent with initial conditions
and we must do it ad-hoc:

ds:=pdsolve({sys}) ;
           {x(t, s) = _F1(s) sin(t) + _F2(s) cos(t),
             y(t, s) = _F1(s) cos(t) - _F2(s) sin(t)}
ds0:=eval(ds,t=0);
              {x(0, s) = _F2(s), y(0, s) = _F1(s)}
subs(ds0, [ic]);
                 [_F2(s) = a(s), _F1(s) = b(s)]
eval(ds,%);
             {x(t, s) = b(s) sin(t) + a(s) cos(t),
               y(t, s) = b(s) cos(t) - a(s) sin(t)}

 

 

For these kinds of maniputations, the vectors, lists, matrices etc were invented.
Example using column vectors:

gamma__w :=10:
z := <2, 1, 0, -1, -4, -7>:
h := <2, 1, 0,  0,  0,  0>:

gamma__w * (h-z);

                   

 

 

 

You have a polynomial systems with many parameters.
The result of solve would be huge and in terms of RootOfs.

Just to make an idea, let's eliminate C4 between the first two equations:

poly:=(numer@lhs)~([eq1,eq2,eq3,eq4]):
resultant(poly[1], poly[2], C4);

You will see that even this first step gives a huge polynomial.

For numerical results, give values to parameters and use fsolve.

 

The residuals are much smaller using LSSolve:

eq:=[seq(lhs(EQQ[i]),i=1..36)]:
Digits:=15:
s:=Optimization:-LSSolve(evalf(eq),iterationlimit=100000);


s := [0.960968723053680e-3, [AA[1] = -.347787136052672, AA[2] = 0.510086502851874e-1, AA[3] = 1.16178642934306, AA[4] = .283780970480705, AA[5] = 0.639701227912076e-1, AA[6] = 0.167348852743906e-2, AA[7] = .970111340350232, AA[8] = 0.123075237465221e-1, AA[9] = -2.56229943109674, AA[10] = .119761186424375, AA[11] = 0.514762128059566e-1, AA[12] = 0.653686059489268e-2, AA[13] = 1.17224901537396, AA[14] = 0.123075188025341e-1, AA[15] = -3.10469172766768, AA[16] = .119761196994863, AA[17] = 0.514762189628706e-1, AA[18] = 0.653685044254999e-2, AA[19] = .970010263800571, AA[20] = 0.122243656427024e-1, AA[21] = -2.56235367446998, AA[22] = .119187577004902, AA[23] = 0.512617524055109e-1, AA[24] = 0.636735006533654e-2, AA[25] = 1.17214793882475, AA[26] = 0.122243606766083e-1, AA[27] = -3.10474597103813, AA[28] = .119187587581231, AA[29] = 0.512617585634991e-1, AA[30] = 0.636733990163037e-2, AA[31] = .516890321334600, AA[32] = 0.510086073253494e-1, AA[33] = -1.15838655820066, AA[34] = .283780982288442, AA[35] = 0.639701289588585e-1, AA[36] = 0.167348506238467e-2]]

Digits:=40:
eval( evalf(eq), s[2]);


Unfortunately LSSolve seems to get stuck for Digits>15.

 

 

Because we already know what has to be proved, it could be done like this:

 

x1 := arcsin(1/2*(2+(2+(2+2^(1/2))^(1/2))^(1/2))^(1/2))

solve(sin(32*x)=expand(sin(32*x1)),allsolutions);

(1/32)*Pi*_Z1

(1)

a:=subs( indets(%)[]=n,%);

(1/32)*Pi*n

(2)

min[index]( [seq(abs(x1 - a), n=1..16)] );

15

(3)

'x1'=eval(a, n=%);

x1 = (15/32)*Pi

(4)

 


 

restart;

sys:= [diff(f(x, y, t, u), u, t)-(diff(f(x, y, t, u), x, y)) = 0,
       diff(f(x, y, t, u), u, u) = 0,
       diff(f(x, y, t, u), u, y) = 0,
       diff(f(x, y, t, u), x, u) = 0,
       diff(f(x, y, t, u), x, x) = 0,
       diff(f(x, y, t, u), y, y, y) = 0];

[diff(diff(f(x, y, t, u), t), u)-(diff(diff(f(x, y, t, u), x), y)) = 0, diff(diff(f(x, y, t, u), u), u) = 0, diff(diff(f(x, y, t, u), u), y) = 0, diff(diff(f(x, y, t, u), u), x) = 0, diff(diff(f(x, y, t, u), x), x) = 0, diff(diff(diff(f(x, y, t, u), y), y), y) = 0]

(1)

S0:=pdsolve(sys);

{f(x, y, t, u) = (_F3(t)*y+_F4(t))*x+(_F3(t)+_C1)*u+(1/2)*_F7(t)*y^2+_F8(t)*y+_F9(t)}

(2)

pdetest(S0,sys);  # bug

[diff(_F3(t), t)-_F3(t), 0, 0, 0, 0, 0]

(3)

S00:=eval(S0, _F3(t)=c*exp(t));

{f(x, y, t, u) = (c*exp(t)*y+_F4(t))*x+(c*exp(t)+_C1)*u+(1/2)*_F7(t)*y^2+_F8(t)*y+_F9(t)}

(4)

pdetest(S00,sys);  # So, S00 is a solution.  But is it the general one?

[0, 0, 0, 0, 0, 0]

(5)

Let us try to answer this question.

 

 

#####################################

 

S1:=pdsolve(sys[2..-1]);  # remove the first PDE

{f(x, y, t, u) = (1/2)*(2*_F5(t)+y*(_F3(t)*y+2*_F4(t)))*x+_F8(t)*u+(1/2)*_F9(t)*y^2+_F10(t)*y+_F11(t)}

(6)

pdetest(S1,sys[2..-1]); #ckeck OK

[0, 0, 0, 0, 0]

(7)

pdetest(S1,sys);

[diff(_F8(t), t)-_F3(t)*y-_F4(t), 0, 0, 0, 0, 0]

(8)

solve(%,_F4(t));

{_F4(t) = -_F3(t)*y+diff(_F8(t), t)}

(9)

S2:=eval(S1,%);

{f(x, y, t, u) = (1/2)*(2*_F5(t)+y*(-_F3(t)*y+2*(diff(_F8(t), t))))*x+_F8(t)*u+(1/2)*_F9(t)*y^2+_F10(t)*y+_F11(t)}

(10)

# S2 should be the general solution of the initial system

pdetest(S2,sys);

[_F3(t)*y, 0, 0, 0, 0, 0]

(11)

# But it is not. Conclusion: the solution S1 of the system S[2..-1] is not the general one.
# So, probably S00 is not general too.
# It should be simple to check this by hand.

###########################################

 


 

Download pde-bug.mw

Edit.  I made a mistake when solving for _F4.
_F3  must be necessarily 0 (due to the presence of y)
So, actually S2 must be

So, indeed, S00 is not general (for example f = t*u+x*y is not contained in S00), but S1 could be after all correct; it would be interesting to check this by hand.
 

 

 

 

subs[eval](F = ( (x,y)->f(x)+g(y) ), F(x,x) );

(Of course [eval] is not mandatory, just to obtain th result evaluated.)

The problem reduces to find/guess  the next term(s) of the sequence:

1, -1/2, 1/12, -1/72, 1/504, -5/18144, 1/27216, -95/19813248

It should be obvious that without extra information this is impossible. Maybe Nostradamus ...

 

Edit. Example:
f := convert( series(cos(x)+x^2, x, 14), polynom):
L := [seq(coeff(f, x, n), n = 0 .. 12)];

    L := [1, 0, 1/2, 0, 1/24, 0, -1/720, 0, 1/40320, 0, -1/3628800, 0, 1/479001600]

rec := gfun:-listtorec(L, u(n));    # not the expected one
     rec := [{(n^3-26*n^2-449*n+538)*u(n+1)+(-707*n^3-1737*n^2+1283*n+393)*u(n+3), u(0) = 1, u(1) = 0, u(2) = 1/2}, ogf]

 


 

ImportMatrix("http://fs3.fex.net/get/245716150875/11071260/data.txt"):
V:=convert(%,Vector);

 

You have generated only a small number of partitions: 1401400 (=nops(P)),  the total number being bell(N-1) = 27644437.

N is too large for a brute force solution (and with your attempt to insert also the permutations, the needed memory is huge).
I think that a very good (or a heuristic) algorithm is needed here.

Optimization:-Maximize((1-y^2)/(x^2),{x^2 + y^2 <= 1, x>=1/2}, x=1/2..1,  y=-1..1  );
            [ 4., [x = 0.5, y = 0.0] ]

 

The answer should be obvious without Maple.
For f(u,v) = 1/(u+v)^2+4*u*v-1,

f(0+, 2) = - 3/4 < 0,  f(1, 1) > 0.
f being continuous in the connected domain (0,oo) x (0, oo),  f must be 0 somewhere.

First 88 89 90 91 92 93 94 Last Page 90 of 120