Robert Israel

6467 Reputation

21 Badges

15 years, 222 days
University of British Columbia
Associate Professor Emeritus
North York, Ontario, Canada

MaplePrimes Activity


These are answers submitted by Robert Israel

Yes, we can... maybe. If you post an example of the sort of problem you're trying to solve, we might be able to help. As Paulina noted in her posting (almost simultaneous to yours), there are two possibilities, the Optimization package and the Global Optimization Toolbox. Important considerations include: 1) How many variables and constraints? 2) Is the problem convex (i.e. maximize a concave function or minimize a convex function, with the feasible region convex)? If not, could there be lots of local optima? 3) Is it a quadratic or least squares problem? 4) How easy is it to find a feasible solution?
For cases where maximize doesn't work, you might try using interval arithmetic (evalr). This gives simple-minded bounds, but that's often good enough for what you need in analysis. For example:
> evalr(eval(2*sin(x)+arctan(x)*cos(x),
     x=INTERVAL(-infinity..infinity)));
INTERVAL(-1/2*Pi-2 .. 1/2*Pi+2)
Since your equations are linear and sqrt(2) is irrational, you can separate the terms that involve sqrt(2) and those that don't. It seems isolve isn't clever enough to do this itself.
> split:= proc(eq) 
    local v;
    v:= lhs(eq)-rhs(eq);
    coeff(v,sqrt(2),1), coeff(v,sqrt(2),0);
  end proc:
  S2:= map(split, S1); 
S2 := {0, kk[3], kk[1], kk[2], kk[4]+kk[5], kk[5]-2}
> isolve(S2);
{kk[3] = 0, kk[1] = 0, kk[2] = 0, kk[4] = -2, kk[5] = 2}
You could do this:
> match(op(2,K) = -op(1,K), x, 's');
  s;
true {b = 1, c = 0, a = 0}
Actually that should be
> evalf(Int(P(x,3),x=0..infinity);
.4583763847 (in Maple 11.01) Numerical evaluation of op(1,g) for large values of z seems to support this. For example:
> evalf(eval(op(1,g), z = 10^6));
.4503764082+0.*I It looks to me like op(1,g) is real for large z. So I think the bug may be in MultiSeries[asympt].
Perhaps you mean something like this. You can use transform in the plottools package to make a 2D plot object into a 3D object (say on a plane in 3D space)
> P1:= plot3d(x^2-y^2,x=-1..1,y=-1..1): # 3D plot
  P2:= plot(x^2 - 1, x =-1..1): # 2D plot
  P3:= plottools[transform]((x,z) -> [x, -1.5, z])(P2):
  plots[display]([P1,P3], axes = box);
The result of isolve is a good starting point. Consider the requirements a[j] >= 0 as a set of constraints on the variables _Z1 and _Z2. Then you could use integer linear programming, available in Optimization[LPSolve] with assume=integer. For example, if S is the result of your isolve, to find the solution that minimizes a[1]:
> Optimization[LPSolve](subs(S,a[1]), map(t -> (rhs(t)>=0), S),
    assume=nonnegint);
[4, [_Z2 = 0, _Z1 = 2]]
OK, I downloaded your worksheet, and it looks just as mystifying as the original message. What are these "⊢" and "⊣"? Are they supposed to be some kind of binary operators? As far as I can tell, Maple doesn't know anything about them. It treats a⊣b as a single symbol.
First of all, the exponential function is exp; e is nothing special. So your definition of f should involve exp(I*w*n), not e^(I*w*n). Second, your equation is not true. To deal with these things in the sense of tempered distributions, we integrate against a smooth rapidly-decreasing test function phi. int(2*Pi*Dirac(w+2*Pi*k)*exp(I*w*n)*phi(n), n=-infinity..infinity) = 2*Pi*Dirac(w+2*Pi*k)*Phi(-w) where Phi is the Fourier transform of phi. int(2*Pi*Dirac(w+2*Pi*k)*Phi(-w),w=-infinity..infinity) = 2*Pi*Phi(2*Pi*k) sum(2*Pi*Phi(2*Pi*k),k=-infinity..infinity) = sum(2*Pi*phi(m),m=-infinity..infinity) (the Poisson summation formula) And this is sum(2*Pi*int(Dirac(n-m)*phi(n),n=-infinity..infinity),m=-infinity .. infinity) so, if I haven't made any mistakes, the correct formula is int(sum(2*Pi*Dirac(w+2*Pi*k)*exp(I*w*n), k = -infinity .. infinity), w = -infinity .. infinity) = sum(2*Pi*Dirac(n-m),m=-infinity..infinity)
If I understand you correctly, you only want to plot the points [t, huur(t)] where t is a positive integer. You can do it this way:
plots[pointplot]([seq]([t,huur(t)],t=0..50));
Use Maple Notation, not 2D Math, for input (Tools, Options..., Display, Input display). No messing around with arrow keys, and you get exactly what you type.
If the calculus book says that, it's actually wrong. The answer is right when |x| < 1. But that doesn't help Maple, which doesn't worry much about convergence of sums. The standard trick (justified by absolute convergence) is interchanging the two summations. AFAIK Maple doesn't have a command for that: it might be a nice addition to the SumTools package. So we do it by hand. Note that the sum is over (i,j) with 1 <= j <= i < infinity, so we take j = 1 to infinity in the outer sum and i = j to infinity in the inner sum.
> sum(sum((-1)^(i+1)/j * x^i, i=j .. infinity), j = 1 .. infinity);
1/(x+1)*ln(x+1)
> `diff/T` := (a,x) -> T(a)*diff(a,x);
Then:
diff(T(g(a)) . h(a), a);


           /        / d      \\                        / d      \
           |T(g(a)) |--- g(a)|| . (h(a)) + (T(g(a))) . |--- h(a)|
           \        \ da     //                        \ da     /
I don't think plottools[transform] would do this, but you could work directly with the Array in the ISOSURFACE plot structure. Something like this:
> plotB := implicitplot3d(f4(alpha,beta,Omega)=0, ...);
# extract the Array
  A:= op([1,1],indets(plotB,specfunc(anything,ISOSURFACE)));
# get its dimensions
  Adims:= ArrayDims(A);
# delta as function of alpha, beta, Omega 
  fdelta:= (alpha,beta,Omega) -> 2*(-beta*Omega/sqrt(2))+Omega/sqrt(2)*
(alpha^2+(3*(abs(beta*Omega))^(3/2)*(IntegralA(1,0,infinity))
/(4*2^(3/4))))*(1/beta);
# Array of delta values
  Delta:= Array(Adims[1..3],
    (i,j,k) -> fdelta(A[i,j,k,1],A[i,j,k,2],A[i,j,k,3]),
     datatype=float[8],order=C_order);
# new Array with Omega, delta, alpha, f4 values
  A1:= Array(Adims, datatype=float[8],order=C_order);
  A1[Adims[1..3],1]:= A[Adims[1..3],3];
  A1[Adims[1..3],2]:= Delta;
  A1[Adims[1..3],3]:= A[Adims[1..3],1];
  A1[Adims[1..3],4]:= A[Adims[1..3],4];
# plot it
  PLOT3D(ISOSURFACE(A1), AXESSTYLE(BOX), 
    AXESLABELS("W","d","a",FONT(SYMBOL,12)));
Your procedure p1 still contains a statement
YP := Array(...)
As I tried to explain before, you should not assign a value to YP. Instead, you should assign values to the entries of YP, which is already an array. Something like
for j from 1 to N-1 do YP[j] := ... end do;
First 132 133 134 135 136 137 138 Page 134 of 138