Functions that are defined by different algebraic expressions over different real intervals are called piecewise, and the corresponding Maple command is piecewise (see help page ?piecewise). Your function could be entered like this (there are several possible variations to this):

f:= piecewise(Or(x < 1, x > 3), 3*(1-x^2)/((x-1)*(x-3)), 3*x);

Another possibility is 

f:= piecewise(And(x >= 1, x <= 3), 3*x, 3*(x+1)/(3-x));

A crude plot can be made by

plot(f, x= 0..5);

Two ways that it can be improved are limiting the vertical range and telling it to look for discontinuities:

plot(f, x= 0..5, -20..20, discont);

You can see from the plot that f is differentiable at x = 1 because the graph has no sharp point there. Symbolic differentiation of piecewise functions also works:

simplify(diff(f,x));

You can achieve more-accurate results for large x by using 

f:= x-> 4*x*ln(1+1/x);

But the round-off error will still win eventually. So, a better solution is to use an asymptotic series:

f2:= unapply(convert(asympt(f(x), x, 10), polynom), x);

It's fairly easy to customize the number of terms in the series to your value of Digits. The series is alternating, so Leibniz's Theorem on the error of truncated alternating series applies: The absolute error is at most the absolute value of the first omitted term.

Entering small matrices is as simple as this (note that the row separators are semicolons, not commas):

A:= <2, -5, 1; 0, 3, -2>;
B:= <0, 2, -1; 7, 1, 3>;
A+B;

First, remove the epsilon[1] and epsilon[2] from the alias command. Having those be functions of x and t will confuse the next step.  They can be put back as functions afterwards if you need. Then do 

Eq2:= mtaylor(lhs(Eq1), [epsilon[1], epsilon[2]], 2) = 0;

The result can be simplified in various ways. But it's linearized with respect to the epsilons.

You are passing to UP4P a procedure f but UP4P is expecting an algebraic expression f. So, you can correct this immediate problem by changing the last argument in the call to UP4P from f to f(x1,x2,x3,x4).

But, why are you doing this? This program is garbage. You should start from scratch, discard the original program, and try to erase it from your mind (because referring back to this program is corrupting your mind regarding good programming practice).

So, what's your goal?

• learning Maple programming?
• learning general computer programming?
• learning numerical analysis (algorithms)?
• learning numerical analysis (programming)?
• just finding minimum and maximum values of functions?

This program is worthless for any of those purposes.

I don't think that you understand what dualaxisplot is supposed to do. Note that the command is passed two plots (and possibly also some "overall" options such as the title in your case---the extra options are irrelevant to my point here). Look at those individual plots. There's no sense in combining the two plots into a dualaxisplot until the individual plots make sense to you. So, do they?

In both of the individual plots, you're asking it to "plot" three numbers for each value of x. What does that mean? 

Also, x is being treated as a continuous real variable in the interval from 1 to 10, not a discrete integer variable. However, this problem is easy to correct if you can answer the questions from the first two paragraphs. 

If you use `f'`, or for better readability `f '`, then the prime will be considered part of the name, and it won't have any mathematical or computational effect. The quotes won't display in output. The characters inside the quotes can be anything, no matter how weird, and the object will be a name.

Like this:

combine(simplify(convert(expr, sincos)))

Yes, it can be easily done with code that is much shorter and more efficient. Let S be that set, or any other structure, to be extracted from. Do

indets(S, 'suffixed'("t", nonnegint) &under String)

The textbook example that you're following starts with two cosine functions of the same period, but different amplitudes and phase shifts. The problem that you're currently trying to solve has a sine and a cosine function of the same period (and different amplitudes and phase shifts). Fortunately, by a trivial identity, sin(x) = cos(x - Pi/2) for any x. Applying that, our two functions are now 

f:= t-> 3*cos(4*t - Pi/2);  g:= t-> 5*cos(4*t - 1);

Applying the textbook's algorithm now, we need the modulus (or abs) and argument (that thing with arctan) of

sc:= 3*exp(-I*Pi/2) + 5*exp(-I*1);

The command polar gives the modulus and argument of any complex-number expression together in one step:

evalf(polar(sc));
             polar(7.697020822, -1.212177299)

And that's the amplitude and phase shift of f + g given as the textbook's answer.

I'm not sure exactly what you mean by "symbol", which is a formally defined type in Maple. Let e be the expression. You could do 

indets(e, symbol)

The problem with that is that it might leave out some things that you want to see. For example, the formal symbols of f(x[n], y) are f, x, y, and n, but the above will only return {y}. The type rules at work here may seem mysterious at first, but it's easy to understand when you know that

1. x[n] is considered a name but not a symbol;
2. symbol is a subset of name;
3. indets doesn't dig down into indices (it does dig all the way to the bottom of expressions, but for some reason indices aren't included), so the n isn't found. 

So, instead, you could do

indets(e, {name, function})

The problem with that is that it may return structures that you think should be further broken down into symbols. For example, if e:= f(x[n], y), it'll return {y, x[n], f(x[n], y)}. You might want that broken down into {f, y, x, n}. That would be fairly easy to do if you tell me that that's what you want.

I don't think that it'd ever be possible to do that in Maple for several reasons. The foremost reason that comes to mind is that in Maple procedures are what some computer scientists call "first-class data structures": They can be manipulated like any other variable. I don't think that that's true in Java. So having a simple variable and a procedure with the same name would be essentially the same as having an integer and a float variable with the same name---and I'll assume that you understand why that couldn't happen in any language.

If by degree you mean the degree with respect to a particular variable, such as the x of your example, then it can be done with simplify by giving an extra argument (or arguments) called side relations, which are polynomial expressions assumed to be equal to 0:

p:= x^6 + 2*x^3 + x^2;
simplify(p, {x^5}); #must use { }
simplify(p, {x^3});

Unlike Christian's Answer, this works even if p is not itself a polynomial but rather is a more-general function of polynomials in x. And unlike Tom's Answer, it also works regardless of whether the polynomials within p are expanded into terms. (However, their Answers, unlike mine, would be fairly easy to modify to handle the multivariate case with degree meaning total degree.)

The side relations facility is much more general than these examples show; see help page ?simplify,siderels.

In your Question, you present the following vector of equations (or expressions implicitly equated to 0) to be solved for the w variables. The equations are linear in those variables.

equ:= < 
    0.85*((phi__cs*beta__s*f__con) . (D__pile*w__1)) + <-5/6*H - 5/8*P>, 
    0.85*((phi__cs*beta__s*f__con) . (D__pile*w__2)) + <5/6*H - 5/8*P>, 
    0.85*((phi__cs*beta__s*f__con) . (D__pile*w__3)) + <1/2*H - 3/8*P>
>;

The form of these expressions is quite awkward. Before proceeding, it's mandatory to correct these 2 things:
1. `.` multiplication must be changed to `*`.
2. Internal vectors must be converted to scalars.

There's another optional correction that I'd like to make, because it makes the output nicer:
3. Convert decimal numbers to exact fractions.

All three corrections can be done automatically, regardless of the number of entries:

eqs:= [evalindets](
    evalindets(
        convert(equ, rational), #Convert 0.85 to 17/20.
        specop(`.`), `*`@op #Correct improper multiplication.
    ),
    rtable, seq #Convert internal vectors to scalars; convert outer vector/matrix to list
);

Your equations are of course trivial to solve individually. However, I suspect that in general you will want to have equations that contain more than one of the w variables each. The code below will also handle that more general case (as long as they're still linear in those variables).

LA:= LinearAlgebra:
W:= indets(eqs, suffixed(w__))[]; #Extract decision variables.
<W> =~ LA:-LinearSolve(LA:-GenerateMatrix(eqs, [W]));

You say that the awkward format of your original equ is due to it being the output of some other program (or equation). But there is no good reason why that program should produce its output in that form. You should change that other program or get its author to change it. 

Constructing matrices from "blocks" (or submatrices) and permuting the rows and columns of matrices is much easier than anything shown in this thread so far, and it can be done using no packages (no LinearAlgebra, no ListTools, no linalg), no loops, and no index variables (such as the i, j, and k that you used).

B1:= 24*<9, -12; -12, 16>; #2x2 block
Z:= Matrix(2,2); #2x2 zero matrix
#3x3 arrangement of 2x2 blocks to produce 6x6 matrix:
K1:= <B1, -B1, Z; -B1, B1, Z; Z, Z, Z>;
K2:= K1[[5.., ..4]$2]; #Apply permutation [5,6,1,2,3,4] to rows and columns.
B3:= <500, 0; 0, 0>;
K3:= <B3, Z, -B3; Z, Z, Z; -B3, Z, B3>;

#This creates a flattened form of your KG:
KG__flat:= <K||(1..3)>;

#If for some strange reason you actually want the non-flat
#KG that you originally had:
KG:= rtable(1..3, 1..1, [K||(1..3)], subtype= Matrix);

Here are some answers to some of your 3 additional questions:

1. How to flatten matrices?
Rhetorical answer: It's much easier to create a flat matrix in the first place (for example, as shown above) than it is to flatten a matrix whose entries are themselves matrices or vectors. So, is there some reason why you can't create them flat in the first place? I'm not saying that flattening is impossible, just that it'd be better to avoid it.

2. How to count the rows or columns of a matrix?
Answer: There are many ways to do it. Here are three ways to count the rows of a matrix M:

[op](1, M)[1];  [upperbound](M)[1];  rhs([rtable_dims](M)[1]);

To count the columns, replace [1] by [2] in any of the above.

3. How to delete rows and columns from a displayed matrix to facilitate the entering of new entries in 2D Input?
Non-answer: I don't know how to do things in 2D Input. However, it's trivial to delete rows or columns, or make them blank, or just re-enter them in 1D input.