restart;
P := <1,1,1> + s*<2,-1,-1> + t*<-1,3,2>;
Q := <1, -3, 5>/sqrt(35);
t_range := 0..1;
s_range := 0..1;
center := eval(P, {t=`+`(op(t_range))/2, s=`+`(op(s_range))/2});
plots:-display([
  plot3d(P, t=t_range, s=s_range, color=blue, style=surface),
  plots:-arrow(center, Q, color=red)
]);

_U1 is a "dummy" variable of integration, because the solution of your differential equation is expressed as a convolution.  You may replace _U1 by any other symbol you like, for example:
restart;
de := diff(x(t),t) = a*x(t) + b(t)*u(t);
dsolve({de, x(0)=alpha}, x(t), method=laplace);
subs(_U1=s, %);
                       

There may be shorter ways of doing this but this is what first came to my mind:

Download worksheet: mw.mw

In the attached worksheet I calculate the fastest path of a point mass sliding without friction on an arbitrary surface z=f(x,y) under the action of gravity which pulls in the negative z direction,. The path's endpoints are prescribed,

For simplicity I assume that the path can be expressed as a function of the x coordinate.  That's  the reason for the "partial" in the title.  The general case requires a fully parametric representation of the path.

Here is an illustration when the surface is z=f(x,y)=x²+y².

In the worksheet there is a formula for the overall travel time, however I don't know how to evaluate it because it is expressed in terms of an integral of an expression involving the solution of dsolve/numeric.  You guys may know how to do that.

brachistochrone-generalized.mw

Note added later

Following Carl Love's suggestion, I have added a section that calculates the travel time.  Here is the extended worksheet.

brachistochrone-generalized-ver2.mw.

The travel time for the illustrative example above is T=4.4019.

I will give a sketch of what needs to be done and leave it to you to code it in Maple.

Consider the rectangle R = { 0 < x < X,  0 < t < T } in the x-t plane.  According to Green's Theorem, for any two functions M(x,t) and L(x,t), the following two integrals are equal:

    P = the integral of (L dx + M dt) on R's boundary;
    Q = the double integral of (dM/dx - dL/dt) over R.

Let's pick L=0, and M=dC/dx. Then we have:
    P = the integral of (dC/dx)dt on R's boundary;
    Q = the double integral of d^2C/dx^2 over R.

1. Let us calculate P
   R's boundary consists of four segments -- two with x constant, and two with t constant.  The segments with t constant do not contribute to the P integral because dt is zero along them.  On the segment with x = 0 we have the boundary condition BC1, that is, dC/dx=0, and therefore that segment does not contribute to the P integral either.  We are left with the segment with x=X.  We conclude that:
       P = the integral of dC/dx (at x=X) from t=0 to t=T.

2. Let us calculate Q
   The PDE in DE1 says that d^C/dx^2 equals dC/dt. (I am taking the diffusion constant to be 1 to simplify my typing, but you may take it anything you wish.)  Then we see that
       Q = the double integral of dC/dt over R
          = int(int(dC/dt, t=0..T), x=0..X)
          = int(C(x,T) - C(x,0), x=0..X).
   This latter integral you can calculate since C(x,0) is the given initial condition, and C(x,T) is calculated by Maple.

We conclude that
    the integral of dC/dx (at x=X) from t=0 to t=T
    =
    int(C(x,T) - C(x,0), x=0..X).

It is not clear what you mean by "Can someone please help me with solving this?".

If you need to test for n=0 or n>0, you already know how to write "if ... then ... end if".  If you need to do something else, then you need to be more specific.

By the way, posting screenshot images is not the best idea.  To help people to answer your questions, please post the actual maple worksheet.  Click on the big fat green arrow in the window where you compose your question in order to upload a worksheet.

It's not clear to me what you are attempting to do (mathematically).  It appears that you have a set of three equations

sys := y = p*x,  x = 100*z/(p*r+r),  z = x*a+y*(1-b) ;

and you want to solve these for x, y, z.  If so, then we do

solve({sys}, {x,y,z});

which produces {x=0, y=0, z=0}.

If you want nonzero solutions, then the coefficients a, b, p, and r cannot be arbitrary.  You will need a relationship like -100*b*p-p*r+100*a+100*p-r = 0.

Enter this in your worksheet:

doit := proc(n::integer)
  local ans;
  if type(n, even) then
    ans := n/2;
  else
    ans := 3*n+1;
  end if;
  return ans;
end proc;

Then doit(4) should evaluate to 2, and doit(5) should evaluate to 16.

These class notes that I wrote some time ago may be of interest to you: inverted-pendulum-notes.pdf

The equation of motion is derived in this Maple worksheet: inverted-pendulum1.mw.
Sorry, I don't have Maple Sim.

Here is an animation of the stabilized inverted pendulum extracted from that worksheet.

This boundary value problem has an exact solution in terms of elementary functions, therefore there is no need for dsolve/numeric.

Details are given in interesting-bvp.mw.

Here is the graph of the solution:

Doing a proc for this is an overkill.  All you need is to multiply the list of the coordinates by Pi/180.  For instance:

L := [36.9,102.04,37.345,101.6104];
Pi/180.*L;

produces [.6440264940, 1.780933969, .6517932092, 1.773436034].

Or even a one-liner:

Pi/180.*[36.9,102.04,37.345,101.6104];

which does the same thing.

Added later:

Or if you want to assign the results to variables, you may do something like:

a,b,c,d := op(Pi/180.*[36.9,102.04,37.34,101.6104]);

After which, a will be 0.6440264940, etc.

I don't know what a "parametric matrix" is but perhaps this will help:

Your reference to xkill leads me to believe that your operating system is linux, and if so, then probably you launch Maple from the commandline.

Maple can freeze when its (undocumented?) JAVAHEAP parameter is too small to handle excessive graphics output. Try increasing it by specifying it on the commandline, as in

maple -x -j 4096

The default JAVAHEAP value is 512. 

The short answer to your question is that yes, if the initial angular velocity is large enough, then the disks will lose contact with the floor.

I have written down the calculations in two-disks-ver2.mw.  You may experiment with the model by varying the parameters.

The diagrams below, extracted from that worksheet, depict the vertical contact force as a function of time in two cases: In the first case the disk rolls somewhat slowly and therefore the contact force remains positive throughout.  In the second case the disk rolls faster and the contact force occasionally becomes negative (which means that the disk detaches from the floor).

Applying the ratio test from calculus we see that the series diverges for all nonzero x, therefore your "unknown" is defined for x=0 only.  There is not much to plot there.