How about something with `applyrule`, such as say,

  map2(applyrule, d^(n::And(posint,satisfies(t->t>=2)))=0, myeq);

I realize that `applyrule` is a command in the Library of interpreted procedures and not a compiled kernel builtin. I had the idea that by "builtin" you meant something else.

acer

The meaning of %1 and related formatting parameters is explained on the error help-page.

The syntax for such parameters in the message string is somewhat similar to the formatting modifiers used in printf, but they are not the same and the explanation will not be obtained simply by looking at the printf help-page.

acer

You can use `rsolve` here, to solve for the general term a(n). And you can use that to produce either a recursive procedure or an expression (in terms of n) for the result.

Note that with this approach you don't first have to adjust (manually) the indices to obtain a(n) = a(0) + a(n-1).

eq := a(n+2) = a(0) + a(n+1):                     

f := rsolve( {subs(a(0)=1,eq), a(0)=1}, a(n), makeproc ):

f(4);                                                    

                                       5

expr:=rsolve( {subs(a(0)=1,eq), a(0)=1}, a(n) );         

                                 expr := n + 1

eval(expr, n=4);                                         

                                       5

It would be nicer if rsolve({eq,a(0)=1},a(n)) worked directly.

acer

You can evaluate your various expressions (psi[1], etc) at the point x=0 without having to assign (or later unassign) to x.

> expr:=A*cos(x)+ R*i*sin(x) = B*exp(-x)+ C*exp(x):

> solve( eval(expr,x=0) );

                           {A = B + C, B = B, C = C}


> expr;

                 A cos(x) + R i sin(x) = B exp(-x) + C exp(x)

Of course that particular example is trivial and doesn't require using `solve`, because,

> eval(expr,x=0);

                                   A = B + C

acer

Simplification with side relations (which Markiyan demonstrates in another Answer) can be an easy way to get such an effect. It's a good answer.

But it can also be noted that there are times and other examples when some of the simplifications done to an expression might not be what is wanted. Consider the case that the input expression is conjugate(z)*z^2*(ln(2*z) and one wishes to replace conjugate(z)*z while also leaving ln(2*z) intact. In other words, the following might not be what is wanted,

simplify( conjugate(z)*z^2*(ln(2*z)), {z*conjugate(z)=x} );

                           (ln(2) + ln(z)) x z

And of course there are related problematic examples where some undesired aspect of simplification is not easily undone after the fact (with `combine` or what have you). In other words, the need for (only) plain substitution is not obviated by fixing the above simplistic example by calling something like `combine` on the result.

So it might be that for some examples one wants just the substitution. And `algsubs` goes awry trying to figure out where and whether z*conjugate(z) is "present". So it's natural to still ask whether `applyrule` could work here.

A difficulty is that `applyrule` does not see z*conjugate(z) as being present within conjugate(z)*z^2 say, and other related situations. An alternative for the rule used might be,

T:=ee->applyrule(z^(m::posint)*conjugate(z)^(n::posint)
                 = x*z^(m-1)*conjugate(z)^(n-1),ee):

And a few rudimentary examples to compare with the simplification with side effects might include,

R:=ee->simplify(ee,{z*conjugate(z)=x}):

f1:=conjugate(z)^2*z^2:
T(f1), R(f1);

                              2   2
                             x , x 

T(conjugate(z)*z^2), R(conjugate(z)*z^2);

                            x z, x z

T(conjugate(z)^5*z^2), R(conjugate(z)^5*z^2);

                              3      3
                           2 _    2 _ 
                          x  z , x  z 

T(conjugate(z)/z^2), R(conjugate(z)/z^2); # no change expected
                             _   _ 
                             z   z 
                             --, --
                              2   2
                             z   z 

But now we might notice that the above does not handle the case where the z*conjugate(z) appears in the denominator. We could try and clean that up crudely, by just repeating the process.

T:=ee->applyrule(z^(m::negint)*conjugate(z)^(n::negint)
                 = 1/(x/(z^(m+1)*conjugate(z)^(n+1))),
                 applyrule(z^(m::posint)*conjugate(z)^(n::posint)
                           = x*z^(m-1)*conjugate(z)^(n-1),ee)):

T(1/(conjugate(z)*z^2)), R(1/(conjugate(z)*z^2));

                             1    1 
                            ---, ---
                            x z  x z

T(1/(conjugate(z)^2*z^3)), R(1/(conjugate(z)^2*z^3));

                            1     1  
                           ----, ----
                            2     2  
                           x  z  x  z

And we can see how the two approaches differ with respect to leaving other terms unchanged (or not), for the example previously mentioned,

expr:=conjugate(z)*z^2*(ln(2*z)):
T(expr), R(expr);

                x z ln(2 z), (ln(2) + ln(z)) x z

It's unfortunate that `applyrule` is not as robust and easy to use as it might be, since it does seem to fulfil at least some needed functionality. (I would not be very surprised if the usage examples I gave above fall down on some related, problematic examples.)

acer

Maple 18 can export 3D surface plots to STL (Stereolithography) file format, which some 3D printers might support directly. (Or you may be able to convert from STL to .x3g format.)

acer

Do you mean like this?

> A:=<<1|2>,<3|4>>;
                                      [1    2]
                                 A := [      ]
                                      [3    4]

> B:=<<2|3>,<4|5>>;
                                      [2    3]
                                 B := [      ]
                                      [4    5]

> A .~ B;          

                                  [ 2     6]
                                  [        ]
                                  [12    20]

> A *~ B;      

                                  [ 2     6]
                                  [        ]
                                  [12    20]

> zip(`*`,A,B);

                                  [ 2     6]
                                  [        ]
                                  [12    20]

acer

It may be that Optimization internal routines are incorrectly identifying the objective expression as having 3 independent variables, where it mistakenly decides that the objective "depends" upon `t` the dummy variable of integration.

You could try something like this (which got me a numeric result in Maple 17.02),

restart:

eta:= 1000:  B:= 5/2:

obj := Int(B*eta^(-B)*t^(B-1)*exp(-(t/eta)^B)*(t-n*h),t=n..(n+1)*h)
       /Int(B*eta^(-B)*t^(B-1)*exp(-(t/eta)^B), t = n .. (n+1)*h):

evalf(subs(n=1,h=1.5,obj));

plot3d(obj,n=0..10,h=1..2);

with(Optimization):

F:=proc(N,H)
  local res;
  res := evalf(subs(n=N,h=H,obj));
  #if not type(res,numeric) then
  #  res:=Float(undefined);
  #  print([N,H]);
  #end if;
  return res;
end proc:

#infolevel[Optimization]:=3:

Minimize(F,0..10,1..2);

Of course you'll have to make sure that the ranges for n and t are correct, and similarly for those ranges of integration (like Carl, I find the range t=n..(n+1)*h a bit strange). Also, you may run into difficulties if it encounters values of n and h for which the numeric integration fails to produce a numeric value for the objective.

I will submit a bug report against the problematic Minimize(obj,n=0..10,h=1..2). I am not sure why it cannot rely upon `depends`, where depends(obj,t) returns false here.

Note that the following also gets a numeric result, by specifying the names of the variables when using the expression form,

Minimize(obj,n=0..10,h=1..2,variables=[n,h]);

acer

Do cylindrical coordinates make your goal easier?

I am just guessing based on your diagram and description, but are you trying to produce something like one of these?

restart:
with(plots):
k:=8:

display(
   spacecurve([1,theta,0],theta=0..2*Pi,
              color=blue,coords=cylindrical),
   spacecurve([1,theta,cos(k*theta)],theta=0..2*Pi,
              color=red,coords=cylindrical,numpoints=1000)
        );

display(
   spacecurve([1,theta,0],theta=0..2*Pi,
              color=blue,coords=cylindrical),
   spacecurve([1+cos(k*theta)/4,theta,0],theta=0..2*Pi,
              color=red,coords=cylindrical,numpoints=1000)
        );

display(
   spacecurve([1,theta,0],theta=0..2*Pi,
              color=blue,coords=cylindrical),
   spacecurve([1+cos(k*theta)/4,theta,sin(k*theta)],theta=0..2*Pi,
              color=red,coords=cylindrical,numpoints=1000)
        );

display(
   spacecurve([1,theta,0],theta=0..2*Pi,
              color=blue,coords=cylindrical),
   spacecurve([sqrt(1+(cos(4*theta)^2-1/2)/2),theta,0],theta=0..2*Pi,
              color=red,coords=cylindrical,numpoints=1000)
        );

display(
   spacecurve([1,theta,0],theta=0..2*Pi,
              color=blue,coords=cylindrical),
   spacecurve([sqrt(1+(cos(4*theta)^2-1/2)/2),theta,sin(8*theta)],theta=0..2*Pi,
              color=red,coords=cylindrical,numpoints=1000)
        );

acer

restart:
with(ScientificConstants):

GetValue(a[0]);
Error, (in ScientificConstants:-GetValue) `a[0]` is not a scientific constant object

GetValue(Constant(a[0]));

                                     -11
                       5.291772081 10   

acer

The SVD may be used on your example. The first r=Rank(C) columns of orthonormal left singular vectors of C may be used, along with the first r singular values, to define matrices A and B.

Here is another diagonalizable 3x3 rank=2 symmetric example, but with eigenvalues of mixed sign.

restart:
with(LinearAlgebra):
C:=Matrix([[1.,1,2],[1,1,2],[2,2,2]],shape=symmetric,datatype=float):

r:=Rank(C);

                                   r := 2

Eigenvalues(C);

                            [-0.828427124746191]
                            [                  ]
                            [                0.]
                            [                  ]
                            [  4.82842712474619]

u,s,vt := SingularValues(C,output=[U,S,Vt]):

map(fnormal, u.u^%T); # orthonormal

                   [1.000000000           0.           0.]
                   [                                     ]
                   [         0.  1.000000000          -0.]
                   [                                     ]
                   [         0.          -0.  1.000000000]

fnormal( Norm(C - u[..,1..r].DiagonalMatrix(s[1..r]).vt[1..r,..]) );

                                     0.

A := u[..,1..r]; # dimension 3x2

                     [-0.500000000000000  -0.500000000000000]
                     [                                      ]
                A := [-0.500000000000000  -0.500000000000000]
                     [                                      ]
                     [-0.707106781186548   0.707106781186548]

B := DiagonalMatrix(s[1..r]).vt[1..r,..]; # dimension 2x3

            [-2.41421356237309  -2.41421356237310   -3.41421356237310]
       B := [                                                        ]
            [0.414213562373095  0.414213562373095  -0.585786437626905]

map(fnormal, B . A );

                        [4.828427125            -0.]
                        [                          ]
                        [        -0.  -0.8284271247]

The repeated eigenvalues {9,9} of your original C form the diagonal Matrix result.

It is also possible to show this on paper, where the difference in sign between columns of u and columns of v account for the differences in sign between the corresponding singular values and eigenvalues.

[update: the above may be a very poor choice of another example, since B.A may not be uniquely determined. comments below, about why, possibly. ugh.]

acer

Is this better?

2014_PDE_project_mod.mw

acer

This can be done without loops.

I wasn't sure whether your various functions names are indexed or not.

f:=M->Matrix(M,M,(i,j)->psi[1,i-1]((2*j-1)/(2*m))):
f(1), f(2), f(3);

f:=M->Matrix(M,M,(i,j)->cat(psi,1,i-1)((2*j-1)/(2*m))):
f(1), f(2), f(3);

acer

f := unapply(a,x);

acer

thingy:=(L1,L2)->remove(member,L1,L2):

thingy([1,2,4,6,2,1,3,6,2],[7,4,2,5,2]);

                        [1, 6, 1, 3, 6]

If the second list is very large and has lots of duplicates not found within the first then it may help to first turn the second into a set. Ie,

thingy:=proc(a1::list, a2::list)
    local s2;
    s2:={op(a2)}; # saves effort, using a uniquified version
    remove(member,a1,s2);
end proc:

thingy([1,2,4,6,2,1,3,6,2],[7,4,2,5,2]);

                        [1, 6, 1, 3, 6]

acer