@vv  A wonderful geometric application of complex numbers.

@Carl Love  In this context, assuming real  is needed not for  solve, but for the correct operation of  Re  and  Im  commands. Therefore, your criticism seems unfounded to me.

Compare:

Re(x+I*y);
Re(x+I*y) assuming real;
 

@Carl Love  Here we do not get anything new. This will not be a solution for any  k , but only for  k=1  or  k=-1 . See:

restart;
sys:={-6*c+(3/2)*c^2-2*b-3*b*c+(3/2)*b^2-3*a*c+(k^2)*b-b+(3/2)*(a^2)+(k^2)*a-a=0,
-2*b-3*b*c+3*(b^2)-6*a*c+2*(k^2)*c-2*c-9*a*b+3*(k^2)*b-3*b+6*(a^2)+4*(k^2)*a-4*a=0,
(3/2)*(b^2)-3*a*c+(k^2)*c-c-9*a*b+3*(k^2)*b-3*b+9*(a^2)+6*(k^2)*a-a=0, 
-3*a*b+(k^2)*b-b+6*(a^2)+4*(k^2)*a-4*a=0, 
(3/2)*(a^2)+(k^2)*a-a=0}:
simplify(eval(sys, {a = (k^2-1)/3, b= 0, c= 0}));
solve(%);

@Sabrina Kamal  Just add this command in the end of your code.

Please upload the worksheet in which this error occurs using the thick green up-arrow in the MaplePrimes's editor.

@Carl Love  Thank you for indicating this inaccuracy.
At first I did not notice the third example at all, but when I noticed, I forgot about this phrase.
By the way, if we are very scrupulous, it is worth noting that, for example,  A = 1 is not always, but if  s<>Pi*k  (k  is an integer) because

eval(A, s=Pi*k) assuming k::integer;
Error, (in assuming) when calling '`one of {cos, eval, sin}`'. Received: 'numeric exception: division by zero'

 

@farah adanan Probably you are right, since in Maple 2018 everything works properly. I do not have Maple 18.

Addition - try to remove this line  local gamma;  as Preben indicated.

@farah adanan Probably you have not fixed everything. Save your worksheet in which the error occurs and load a link to it here using the bold green up-arrow in the mapleprime's editor. I'll check.

@Markiyan Hirnyk  Obviously the reason is in  rational  option. Examples show that in a solidly filled region the function takes complex values.
Look at these simple examples:

plots:-implicitplot(sqrt(-x^2-y^2), x=0..1,y=0..1);
plots:-implicitplot(sqrt(-x^2-y^2), x=0..1,y=0..1, rational);
plots:-implicitplot(sqrt(-x^2-y^2), x=0..1,y=0..1, gridrefine=2, rational);
 

@nm  The reason is probably in this:

whattype(1/sqrt(x));
op(1/sqrt(x));
                                                

A workaround:

expr:=1/sqrt(x)+sqrt(x);    
subs([x^(1/2)=abs(x)^(1/2), x^(-1/2)=(abs(x))^(-1/2)], expr);
                                           

 

@digerdiga  Yes. Here is an example:

simplify(sqrt((R-sqrt(R^2+z^2))^2));
simplify(sqrt((R-sqrt(R^2+z^2))^2))  assuming R-sqrt(R^2+z^2)<=0;
                              

@davtrs  This works. I added different colors to your lines.

Download Plot_problem_new.mw

@Carl Love I just thought  D(x)  was just a name, for example:

D(x):=5:
D(x);
                                     5

@vv 

In the original code I replaced  0.5  by  1/2 :
is(D(x)=DD(x));
                                        true

@whtitefall  Due to the parity of the integrand function, the answer to your question can be given without any computation  x=0 . Unfortunately, it seems that you have not understood my answer at all and again you have a lot of mistakes in your code.