@mehran rajabi The usual way for searching for roots is to plot the function in order to get a rough idea about where the roots are.  Then call fsolve with appropriate ranges to pinpoint them.  In your case we do:

plot(eq, x = 50 .. 200, discont = true, view = -100 .. 100)
fsolve(eq, x = 60 .. 70);
                          68.93033112
fsolve(eq, x = 70 .. 120);
                          106.9740567
fsolve(eq, x = 120 .. 170);
                          151.0624355
fsolve(eq, x = 170 .. 210);
                          198.3826229

@Adjal yassine I am quite familiar with the finite element method and its uses for solving differential equation, especially those that relate to structures. I am puzzled, however, by your statement that you are interested in the stiffness matrix, but not in the solution.  What is the use of a stiffness matrix other than for computing a solution?

Perhaps you have good reasons for writing your programs but those reasons do not come through in your posts.  If you wish to attract people's attentions, you should be more explicit in explaining your goals.

@acer That's excellent.  I have your setsize in my ~/.mapleinit now.  It works even with animations. Thanks!

The presentation in your worksheet makes things more difficult than is necessary.  Consider the following.

1. You have defined the variables y, B, R, R2, r as functions.  To make life easier, I suggest making all of them into expressions.  For instance, replace

y := (x,k,xi,B) -> whatever

with

y := whatever

2. In the definitions of R and R2 replace B(k) with B.

3.  Is the y that enters the definition of r the same y that is defined on the first line in the worksheet?  I don't think so.  And if that is indeed a different y, then you should name it something different.

4. In order to plot any of these functions, you need to specify the value of k.

@CEFOG Look at y = a*x + b.  Can you solve for y/a so as to get y/a = (something that has no a in it)?  Do it by hand, not Maple. If you show me how you do that by hand, I will show you how to do it in Maple.

@Joa The symbolic solution that you have obtained is too unwieldy to be useful.  At least I cannot make it do what you want.  But the problem that you wish to solve is a well-explored subject and much has been written about it.

Specifically, we have a system of differential equations which involve a certain number of parameters.  We wish to determine the parameters so that the solutions fit a prescribed set of data.  Doing that is not very difficult but requires some familiarity with common techniques in numerical computing.  If you have a serious need (such as graduate research) for solving the problem that you have stated, I suggest that you have a look at Notes on Adjoint Methods.  Section 6 addresses exactly what you need.

@Earl Your calculations are correct, however I would organize the worksheet somewhat differently.  You may continue doing it your way, or do it as shown below, whichever makes better sense to you.

Download Slinky2-alt.mw

@Earl  Page 212 of that book reduces the solution of the problem to solving a system of differential equations for the unknown functions r(m*) and F(m*)  [ equations (3) and (4) ]  subject to the boundary conditions in (2).

The rest of that page muddles things up by attempting to establish an analogy with harmonic oscillations.  The analogy is quite confusing since there is no time-dependence in this problem—we are looking at the steady-state rotation of a stretched spring (unvarying length) about an axis.  There is not much to animate here.

Here I will show how to solve the equations noted above in Maple.  The note at the top of page 213 does the same thing, therefore the confusion of page 212 was not necessary at all.
 

Download mw.mw

@imparter I know that alpha is a parameter that you want to vary.  The problem is, the code that you have shown has a few dozen of ohter symbols whose meanings are known only to you and the author of the book that you are reading.  I can only guess what those symbols mean.  In fact, I can only guess what the code is attempting to solve because you have not even provided me with the mathematical statement of the problem.

That said, I have made a best effort to guess what you have in mind and have produced this worksheet:  fixed-alpha.m which plots the solution curve at t=0.2 as you have asked, but for alpha=1 only.

To plot the curve for other values of alpha you have two choices:

Choice #1. Copy and paste the entire code (but not the "restart" line) in the worksheet (thus doubling the number of lines in the worksheet).  Set the value of alpha in the second half, execute the entire worksheet.  Then you will have two graphs, the original one (with alpha=1) and the second one with whatever value of alpha you want. 

Copy/paste some more to produce additional graphs with other choices of alpha.  Then you may produce a composite graph by applyting plots:-display.

Choice #2. Instead of copying/pasting, you may wrap the entire code in a proc which takes alpha as a parameter.  Then you may execute that proc any number of times with the desired values of alpha.  This second approach is the right way of doing it but it needs some elementary knowledge of Maple programming.  I cannot tell how much of Maple programming you know.  If you don't, then go with Choice #1.

@imparter The plot the curves based on the code that you have shown requires reading the book to understand what the various variables mean. But I don't have the book. Sorry.

@imparter In my previous message I showed how to plot the solution at time t=0.2 for alpha = 0.1, 3.0 and 5.0. Those are the curves drawn in red, green, and blue.  Isn't that what you are asking?

@AhmedRahby The first step in changing the order of integration, is to draw an accurate picture of the domain of integration as I did in my answer to your previous question.  Did you see that?

It can't be done without pictures.

@Carl Love Thanks for the complement. I, too, teach multiple integrals that way.  Some students see the idea right-away and have an easy time.  For some others it never clicks because they have difficulty in connecting analysis with the underlying geometry.

@AhmedRahby  Do you realize that a multiple integral expresses integration over a certain domain?  The integration limits specify the domain.  Can you sketch the integration domain for your integral?  There is no point in playing around with integration if you don't see that connection.

You have terminated some of the crucial parts of the code with colons and therefore deprived yourself of the chance to visually inspect what Maple is doing.  Don't do that!  Terminate commands with semicolon and verify that Maple echos what you expect.

If you change colons to semicolons, you will see that pde1 contains the undefined variable ρA and ibc1 contains the undefined variable ℓ. Correct those!

After those corrections your problem is well-posed and should have a solution.  Unfortunately Maple's pdsolve is unable to find it.  The fact that we have v(t) rahter than v(x,t) is not the cause of the failure (see the worksheet below.)  The problem arises due to the terms w(0,t) and w(ℓ, t) that occur in pde2 which to Maple looks like a boundary condition rather than a PDE, and therefore it complains about insufficient number of equations.  I don't see a workaround.  If you really need a solution, you may write your own finite difference solver for it.

The much simpler system in the worksheet below suffers from the same issue.

Download mw.mw