Complete code - where I just guessed a value for C_A0 (so change it to anything you want) - would be

restart;
sol1:= dsolve
           (  [ diff(X(W), W) = (0.536000000000000e-3*(1-X(W)))*(1+X(W)),
                X(0) = 0
              ],
              numeric
           ):
C_A0:= 5.0:
C_A:= C_A0*(1-X(W))*(1+X(W)):
plots:-odeplot( sol1,
                         [  [W, X(W), color=red],
                            [W, C_A, color=blue]
                        ],
                        W=0..10000
                     );

The only stl file I have access to is the one in Maple's data/example directory. It is relatively simple to extract the polygons from this file, and then (collectively) animate them.

See the attached

animstl.mw

I'm pretty unfamiliar with the Physics package but I manage to get an output from your last execution group, after fixing (what seemed to me) a few syntax-type things. Try

Trace(-LeviCivita[mu, nu, rho, sigma]*Ucn(X).(((1/6)*e*A[nu](X)+M5*KroneckerDelta[nu, 0])*~Matrix(2, 2, [[1, 0], [0, 1]])).Ucn(X).d_[rho](Ucn(X)).Ucn(X).d_[sigma](Ucn(X)));

Make sure that the directory exists: in youor filename path, you have

"C:/Users/test/Desktop/txtop1.txt"

Is there really a User called 'test' on your system. Even if this directory exists, I'm not convinced about the syntax of youu print statement.

See the attached, where I have added an execution group to perform the function you require. This works on my machine, but you will have to change the path of fname1 to a directory which exists on your machine

pdePrint.mw

whihc produces the result file

txtop1.txt

odeplot() allows you to plot more or less any function of the quantities returned by the dsolve() command. In your case it is relatively simple to isolate an expression for diff(sigma(t),t) from the original ODE, given as 'nde'. Hence one can use odeplot() to plot this expression.

Data can be retrieved from plots using the plottools:-getdata() commands. The resulting matrices can then be exported (in a variety of formats) using the ExportMatrix() command

See the attached

odeProb.mw

Try

seq(is(tm[i] < (1/6)*Diameter[rør][i]), i = 1 .. 17);

or

map(is, [seq(tm[i] < (1/6)*Diameter[rør][i], i = 1 .. 17)]);

but this (still will not return a solution for x,y! Consider

A-B+x^y/y=0;
x^y-y=0

Now the second equation requires x^y=y. Make this substituion in the first equation to get A-B+y/y=0, ie A-B+1=0. In other words your equations are conssitent for any values of x ,y provided A=B-1

If you actually want 3D plots, then you need to use DEplot3d().

Some possibilities are shown in the attached

dePlots.mw

rlLopez response, ie

restart;
P := x + y + a*x^2 + b*x^3 + c*x*y + d*x^2*y + e*x*y^2 + y^2:
mtaylor(P,[x,y],3) - mtaylor(P,[x,y],2);

returns

a*x^2+c*x*y+y^2

when it should return

a*x^2+c*x*y+y^2+x+y

Note the last two entries in the above!

Something like

restart;
P := x + y + a*x^2 + b*x^3 + c*x*y + d*x^2*y + e*x*y^2 + y^2:
select( i->is(degree(i,x)+degree(i,y)<3), P);

returns

a*x^2+c*x*y+y^2+x+y

and is probably a better bet.

select/remove commands are very useful for this kind of exercise. For example one could achieve the same result by removing all terms with degree >=3 (rather than keeping all terms with degree<3) using

restart;
P := x + y + a*x^2 + b*x^3 + c*x*y + d*x^2*y + e*x*y^2 + y^2:
remove( i->is(degree(i,x)+degree(i,y)>=3), P);

1. You could try reading the help page for  'for' by entering ?for at the Maple prompt
2. You could try reading the help page for  'if' by entering ?if at the Maple prompt
3. Given an inability to read help pages, you might consider the following
   1. Maple allows the "block" assignment of contiguous entries in an indexable entity - so actually you do not need any 'for' loops at all.
   2. It seems as if you indexable entity 'w' has zero-based indices (eg you refer to w[0,0]) nut all of your 'for' loops start at 1. This may (or may not) be deliberate

Bearing in mind the above (particularly item (3) ), you might want to consider the attached

nonLoops.mw

return equations?

return assignments?

I don't get it!

Maybe something in the following will help?

returns.mw

See the attached worksheet where I have used only the standard 'methods' returned by pdsolve().

I have been unnecessarily "verbose" in the attached, in the hope of clarity. This is probably(?) unnecessary but since I am not quite sure exactly what the issue is here - I can't be sure

pdeProb.mw

but you are going to have to work with me here. I can see all sorts of places whre this can be "tidied-up", but since I have no idea about your "ultimate objective" I'm going t have to go through this pretty much line by line.

I am not going to get all the way to the end of your code on a "first pass" - we start with the first few lines, and (hopefully) work on from there!

Your first entry

Stiffness:=Matrix(6,6):
Stiffness:=Matrix(6, 6, {(1, 1) = 20.4, (1, 2) = .537, (1, 3) = 1.01, (1, 4) = 0., (1, 5) = 0., (1, 6) = -0., (2, 1) = .537, (2, 2) = 1.59, (2, 3) = .73, (2, 4) = 0., (2, 5) = 0., (2, 6) = -0., (3, 1) = 1.01, (3, 2) = .73, (3, 3) = 2.86, (3, 4) = 0., (3, 5) = 0., (3, 6) = -0., (4, 1) = 0., (4, 2) = 0., (4, 3) = 0., (4, 4) = 1.5, (4, 5) = 0., (4, 6) = -0., (5, 1) = 0., (5, 2) = 0., (5, 3) = 0., (5, 4) = 0., (5, 5) = 2., (5, 6) = -0., (6, 1) = 0., (6, 2) = 0., (6, 3) = 0., (6, 4) = 0., (6, 5) = 0., (6, 6) = .8});

You don't have to "index" each entry separately. There are many ways to shorten Maple's matrix input: all of the shorter versions require about the same amount of typing, so I am just going to show one. The above is entirely equivalent to

Stiffness:=Matrix( 6,
                              6,
                              [  [  20.4,        .537,  1.01,  0.,   0.,  -0.   ],
                                [       .537,  1.59,      .73,  0.,   0.,  -0.   ],
                                [     1.01,      .73,    2.86,  0.,   0.,  -0.   ],
                                [     0.,       0.,         0.,     1.5,  0., -0   .],
                                [     0.,       0.,         0.,     0.,    2., -0.   ],
                                [     0.,       0.,         0.,     0.,    0.,    .8 ]
                              ]
                           );

You may thing that the above code contains an excessive amount of indentation/spacing/whatever, but I write stuff this way because I find it much easier to read. If you want to write the above as

Stiffness:=Matrix( 6,6,[ [ 20.4,.537, 1.01, 0.,  0., -0.],[.537, 1.59,.73, 0., 0., -0.],[1.01, .73, 2.86, 0., 0., -0.],[ 0., 0., 0., 1.5, 0., -0.],[  0., 0., 0., 0., 2., -0.],[0., 0., 0., 0., 0., .8 ]] );

then that is absolutely fine with me! Either of the above solutions is waaaaay shorter than typing your original.

Next point is rotation matrices: you really, really do not have to type these explicitly, because given an axis and an angle, Maple will generate them for you ( you can check out the help for these by entering ?rotationmatrix at the Maple prompt), but basically all you need to know is the rotation angle, and the axis of rotation, so consider the following examples

with(Student[LinearAlgebra]):
R1a:=RotationMatrix(Pi/4, <1,0,0>); # rotation  Pi/4 about +ve x-axis
R1b:=RotationMatrix(Pi/4, <-1,0,0>); # rotation  Pi/4 about -ve x-axis
R1b:=RotationMatrix(-Pi/4, <1,0,0>); # rotation  -Pi/4 about +ve x-axis
R2:=RotationMatrix(Pi/4, <0,1,0>); # rotation  Pi/4 about +ve y-axis
R3:=RotationMatrix(Pi/4, <0,0,1>); # rotation  Pi/4 about +ve z-axis

So don't write ytour rotation matrices explicitly - so long as you know the angle, and the axis of rotation, then Maple will generate the matrices for you. Doesn't matter what the angel is or what the axis is, Maple will still produce the matrix, so for the following totally random angle, totally random axis, Maple will still generate the relevant matrix

Rx:=RotationMatrix(10003.754, <23.7,-15.8,17.1>);

Now you don't seem to want these rotation matrices - rather you want their transposes (no idea why!?), but then just "nest" the commands, in order to save typing, as in

with(Student[LinearAlgebra]):
T1:=Transpose(RotationMatrix(Pi/4, <1,0,0>)); # transpose of rotation  Pi/4 about +ve x-axis
T2:=Transpose(RotationMatrix(Pi/4, <0,1,0>)); # transpose of rotation  Pi/4 about +ve y-axis
R3:=Transpose(RotationMatrix(Pi/4, <0,0,1>)); # transpose of rotation  Pi/4 about +ve z-axis

Now I'd really like to shorten up the definition of your matrix R. I'm assuming that there is some kind of function/pattern, because if I could work out the general rule for how R[i,j] depended on the entries of S5, then I could make this definition much neater. But I've stared at this at your definition for a while and I can't "see" the pattern/function - so if you want this in a "short/clever" way, you are going to hav to tell me what the "pattern/function" actually is!

The attached solves three different ODE systems, and plots output

First system uses a[1], b[1], c[1], d[1].
Second system uses a[2], b[2], c[2], d[2].
Third system uses a[3], b[3], c[3], d[3].

ode.mw

Your code has

c:=[0,1];

therefore 'c' is a list. This is used with

element(s,c);

but arguments to the procedure 'element' are checked with

......proc(s::list(set), c::array(1, nonnegint))

so the second argument must be an array. You are passing a list so the parameter check will fail.

The following will work

s:=[{alpha,beta,gamma},{x,y}];
c:=array([0,1]);
element := proc(s::list(set), c::array(1..nonnegint))
    local i, j;
    for i to nops(s) do
        c[i] := c[i] + 1;
        if c[i] <= nops(s[i]) then
            return [seq(s[j][c[j]], j=1 .. nops(s))];
        end if;
        c[i] := 1;
    end do;
    c[1]:=0;
    FAIL;
end proc;
element(s,c);

but you should be aware that array() has been deprecated for a loonnngggg time. It may be used in the book to which you refer, but if you are using a reasonably up-to-date Maple version then such use *may* have issues