@digerdiga That I cannot reproduce the error for ftoc and ftocms for f probably just means that I made a mistake when I wrote in my answer that " For f the "successful" methods are ftoc and ftocms as for fc, but in addition the correct meijerg.".
Whas I take back is the claim that ftoc and ftocms gives the same results for f and fc.

The results from ftoc and ftocms for fc are equal and wrong.
You use Maple 2015 I see. Also there the results for f are correct for all the three metods ftoc, ftocms, and meijerg as is shown by the result of
 

int(f,t=0..infinity,method=_RETURNVERBOSE) assuming s>1;

which is

@digerdiga FTOC: Suppose F is an antiderivative of the continuous function f on an interval I=[a,b], i.e. F'(x) = f(x) for all x  in I. Then int(f(x),x=a..b) = F(b)-F(a).
Basically I believe that is the method used: Find F and evaluate at the ends of the interval.
But the function F actually found by Maple may not be an antiderivative on all of the given interval (here [a,b] ), but only piecewise so. Thus Maple (or you) will have to change this F into a continous F1 as illustrated in this example:

f:=1/(cos(x)+2);
F:=int(f,x);
simplify(diff(F,x)); # OK generically
plot(F,x=-3*Pi..3*Pi,discont=true); # discontinuous
F1:=int(f,x=0..xx) assuming xx>Pi;
plot(F1,xx=0..3*Pi);

Another problem is that F may not be defined at one or both endpoints of the interval. Thus eval(F, x=a) must be replaced by limit(F, x=a, right) and similarly at the other end point.
An example of that is:

f:=1+ln(x);
int(f,x=0..2);
F:=int(f,x);
eval(F,x=0); #Error
limit(F,x=0,right);# Fine

@Carl Love Surely you are aware of the following special behavior of evalf, but I shall use the opportunity to point it out nevertheless.
evalf accepts a sequence if not actually typed in as arguments one at a time.

seq(Pi^k,k=1..7);
evalf(%); #OK
evalf(seq(Pi^k,k=1..7)); #OK
evalf(seq(Pi^k,k=1..0,-1)); #OK 
seq(Pi^k,k=1..0,-1);
evalf(%); #OK
evalf(Pi,1); # the old alternative to evalf[1](Pi)
evalf(1,2,3,4); # Error

@Carl Love The wrong integral has fc as integrand. So you ought to add the corresponding lines for fc. But the result is indeed correct:
 

fc:=evalc(f); # Your f
int(fc, t= 0..infinity) assuming s > 1, r > 0:
simplify(eval(%, r= 1/2));

@digerdiga The problem here is that the indefinite integral of fc is wrong. That surprised me.
 

restart;
f := ((1/2-I*t)^(-s)-(1/2+I*t)^(-s))/(2*I);
fc:=evalc(f);
## Checking that fc = f for s=2 on t=0..5:
plot(eval(f-fc,s=2),t=0..5); #OK zero
##
Fc:=int(fc,t) assuming s>1, t>0;
Fcdiff:=diff(Fc,t);
## Now Fcdiff and fc ought to be equal:
plot(eval(Fcdiff-fc,s=2),t=0..5); # Not zero
## ftoc and ftocms uses the Fundamental Theorem of Calculus.
## Thus the erroneous result is due to Fc being wrong:
limit(Fc,t=infinity) - limit(Fc,t=0,right) assuming s>1;
simplify(%);

The following session where printlevel is just the default 1 also makes the result 0:
 

restart; 
interface(typesetting=extended); 
debug(`shake/shake`);
evalf(frac(Pi^20)); # 0
printlevel; # 1

Also here the result is 23 if typesetting = standard instead.

Since the second session shows an error from Typesetting:

Error, (in j) invalid subscript selector

<-- ERROR in \`shake/shake\` (now in \`evalr/shake\`) = unexpected result from Typesetting}

I tried suppressing output as in
 

restart;
interface(typesetting=extended);
printlevel:=40:
evalf(frac(Pi^20)):
%;

The error didn't show up and the result was 23.

@farah adanan I tried your worksheet as given in the link. I used Maple 18.1. I found no problem.
I also checked Maple 15, 16, and 17.
The only problem appeared in Maple 15. It didn't like the local declaration at the start. That feature wasn't introduced yet.
The way you and Kitonum are using parameters via e.g. eval(eq1, fixedparameter1), the local declaration can be left out anyway since gamma is being replaced by the value 1 given in fixedparameter1.

Very nice.
Two comments on the worksheet itself:
1. You are apparently setting kernelopts(opaquemodules=false). Otherwise `pdsolve/BC`:-EvaluateBCAtSolution wouldn't work.

2. You have a line which in my downloaded version looks like simplify( Equality - Equality), which to the uninitiated obviously is zero and should be without simplify. Apparently some labelling is going on which I don't see, probably because I don't use labels ever.

@mmcdara Actually the result of dsolve(sys) doesn't look bad if the RootOfs are replaced by aliases (or just by names).
It seems that Im can decide if the roots areal.
 

restart;
sys := {
diff(tau__1(t),t)=phi__1(t),
diff(tau__2(t),t)=phi__2(t),
diff(phi__1(t),t)+2*phi__1(t)+3*phi__2(t)+4*tau__1(t)+5*tau__2(t)=6,
diff(phi__2(t),t)+8*phi__1(t)+7*phi__2(t)+9*tau__1(t)+10*tau__2(t)=11
};
sol:=dsolve(sys);
RO:=indets(sol,specfunc(RootOf));
Im~(RO); # {0} so all roots are real
Re~(RO); # Just returns the input consistent with all roots being real

I don't use alias very often, but here it might be useful:
 

alias(seq(alpha[i]=RO[i],i=1..4));
sol; # Short
evalf(sol);

Be aware of the following point in the help for alias:

"Because aliases are resolved at the time of parsing the original input, they substitute literally, without regard to values of variables and expressions.  For example, alias(a[1]=exp(1)) followed by evalf(a[1]) will replace a[1] with exp(1) to give the result 2.718281828, but evalf(a[i]) will not substitute a[i] for its alias even when i=1.  This is because a[i] does not literally match a[1]."

This implies that the following attempt to find floating point values for the roots won't work:
 

seq(alpha[i],i=1..4);
evalf(%); # no floats, just names
## The simple way works
RO;
evalf(%);

Ordinary names could be used instead of aliases as done here:
 

alias(seq(alpha[i]=alpha[i],i=1..4)); #Removing the aliases first
RO; 
#sol; # Now big
sol2:=subs(seq(RO[i]=alpha[i],i=1..4),sol); 

sol2 is as short as sol with alias, but of course evalf won't do anything about the roots in sol2.

A last comment: I tried evalc  in order to get an expression that also on the surface looks real:
 

r1:=allvalues(RootOf(_Z^4+9*_Z^3+4*_Z^2-19*_Z-5, index = 1));
r1a:=evalc(r1);

I succeeded in some attempts (after restarting) and other times it took too long for my patience.
By using indets(r1a,function) I found these:
{arctan(3/3268*sqrt(10275765)), cos((1/3)*arctan(3/3268*sqrt(10275765))), cos(2/3*arctan(3/3268*sqrt(10275765))), sin((1/3)*arctan(3/3268*sqrt(10275765))), sin(2/3*arctan(3/3268*sqrt(10275765)))}
All very good, but indets(r1a,radical) brought some square roots of not oviously positive quantities.

@annarita To get that error you must be using 2D input.
In 2D input a space is interpreted as multiplication. My code had a space between RGB and ( R() , R(),R() ).
While that is fine in 1D input (aka Maple input) it isn't in 2D input.
Just remove the space so that you have
 

 RGB( R() , R(), R() )

I tested your code in Maple 12 and Maple 15 besides in Maple 2018 and didn't run into any problems. I don't have Maple 13 on this machine.
The 3 results agreed after simplification.

@Carl Love odetest works fine on the implicitly given solutions regardless of the names or values of the constant.
On the explicit solutions having the RootOf odetest converts to an implicit form, which is seen by setting infolevel[odetest] to 2 or higher:
 

restart; 
ode := diff(y(x), x) = x*ln(y(x)):
sol:=dsolve(ode, implicit);
odetest(sol,ode); # 0
odetest(subs(_C1=C1,sol),ode); # 0
odetest(subs(_C1=1,sol),ode); # 0
explicit_sol:=op(solve(sol,{y(x)}));
infolevel[odetest]:=2:
odetest(explicit_sol,ode); # 0
odetest(subs(_C1=C1,explicit_sol),ode); # problem
odetest(subs(_C1=1,explicit_sol),ode); # problem

The procedure used for the job is `odetest/implicit`.

@pppc I used Q as a constant and should have mentioned it above.
The method still works for Q = (1/2)*q*(x-50)^2-1250*q, however, but you need to give q a value before plotting, e.g. q=1 as in

plot(eval(U, {C__1 = 0.12e9, C__2 = 0.2e10, C__3 = 0.2e9, C__4 = 0, C__5 = 0.3e8, q = 1}), x = 0 .. 100, size = [2000, default]);

@guto_vasques You simply forgot to write range = ... . Surely you meant to write:
 

p2b := dsolve(sys2, type = numeric, abserr = 1.*10^(-12), relerr = 1.*10^(-12), range=te .. tr);