@Rouben Rostamian  

I don't understand the mathematics here, so this comment may be irrelevant: Before you get too far involved with this problem in a Mapleish way, remember that Maple currently will only numerically solve PDEs with exactly two independent variables. There's no such hard limit for symbolic solutions.

You were using Maple's 2D input and you typed the epsilon in the first line too close to the following left parenthesis, so it was interpretted as if epsilon were a function rather than as a coefficient. If you use Maple's 1D input, this and numerous other problems will be avoided. In the second line, you will get better results using eval rather than subs. The eval forces the evaluation of the derivatives. So, change these two lines to

eq:= diff(x(t),t$2) - epsilon*(alpha*x(t)^4 - beta)*diff(x(t),t) - x(t) + x(t)^3;
eq2:= eval(eq, x(t)= 1 + epsilon*x[1](t) + epsilon^2*x[2](t));

@Christopher2222 How does one hide a code edit region? I see the option to collapse it, but that doesn't completely hide it.

@MDD 

If you need for your procedure to be able to detect empty Matrices, then I recommend looking at the ColumnDimension rather than the RowDimension, because, like I said, any reasonable attempt to create an empty Matrix out of square brackets leads to a Matrix with one row and zero columns. (<<[][]>> creates a Matrix with zero rows and one column.)

If you need for your procedure to create empty Matrices, then use Matrix() or, if speed is critical, rtable(1..0, 1..0, subtype= Matrix).

If you post your procedure, I can make further suggestions.

@Haley_VSRC 

The corresponding mathematical reason is that b must be less than or equal to the minimal value of t in order to avoid nonreal values. The graphical result is that the omitted points are very close to the regression curve horizontally.

@Markiyan Hirnyk 

No, the {} indicates that the set of elements is the empty set, which is different from the empty set being an element.

@MDD 

NonzeroRowsAndColumns:= (M::Matrix)-> [
     add(`if`(ArrayTools:-IsZero(M[j,..]), 0, 1), j=1..op([1,1], M)),
     add(`if`(ArrayTools:-IsZero(M[..,j]), 0, 1), j=1..op([1,2], M))
]:

@Markiyan Hirnyk 

When I said that my fit had lower residuals, I was using the fit that treated the first four points as outliers (my Sol_Carl5), and I compared it to doing the equivalent with DirectSearch:-DataFit, using all nine fitmethods and taking the minimal residuals. Done this way, my method has less than half the residuals.

But, if I use all four methods crossed with all nine fitmethods (for a total of 36 calls to DataFit), I can make a small improvement over my method. Here's the worksheet:

Download RecenteredPowerFit_vs_DirectSeach.mw

@Markiyan Hirnyk I'm not trying to blame anybody! I was just trying to tell you that both issues were caused by the same mistake. Geez, I can't even say that you're right without you finding some way to criticize!

@capea 

You didn't execute the expandoff command! That is nescessary to prevent simplification of residues. Also, don't forget to restart before the whole process.

@Markiyan Hirnyk You're absolutely right. I made a mistake. I will correct the Answer now.

Regarding the expand portion, you're experiencing the same mistake. After I make the correction to the Answer, try expand(Psi(3)) before and after the expandoff commands.

@MDD Try using the Maple help: Enter ?Psi at the Maple prompt.

Psi(z) = diff(ln(GAMMA(z)), z)

@MDD You simply cut and paste and execute that code that Roman gave, then 

with(PolynomialIdeals):
EliminationIdeal(<x>,{a,b,c});

@Markiyan Hirnyk 

By what criterion are you saying that the result produced by DirectSearch:-DataFit is better than that produced by my RecenteredPowerFit (which calls Statistics:-PowerFit)? You can't directly compare the residuals as reported by the various fitmethods because they are each computed by their own formula. If you compare the sum of the squares of the residuals, my fit is several times better than yours. It is also clear from the plots that my fit is better than yours. (My fit is the red curve in my plot; ignore the blue curve.)

Any result with c > 1 is completely ridiculous, not just "not so good". Try plotting them. Of the nine fitmethods, seven produce such results for me (on the first run).

Furthermore, the results of DirectSearch:-DataFit are not entirely reproducible. Rerunning the exact same call (of course, without restart) often gives a completely different result, often changing a decent result into a completely ridiculuous one.

The vote up may be from the OP and may be due to the fact that I worked back and forth with her over three days and that I produced code that she was able to execute and with which she was able to get what she considered to be "perfect" results.

@Jerome BENOIT I also don't care what the shown material is called. I just didn't want other readers to get the impression that it was ever possible to see the complete code of a module---not by eval or by any other means. Your usage of "body" was likely to give that false impression.