@Hong-Yee Chiu

Please make sure you know the rules: https://wolfr.am/READ-1ST

The rules explain how to format your code properly. If you do not format code, it may become corrupted and useless to other members. Please EDIT your post and make sure code blocks start on a new paragraph and look framed and colored like this.

int = Integrate[1/(x^3 - 1), x];
Map[Framed, int, Infinity]

Thanks for your suggestion. It has problems as follows. (1) The original soln1 gives the following results:

soln1[a_, b_, c_, \[Xi]_, \[Lambda]_, \[Phi]_, rn_] := 
  NSolve[\[Rho]^2 + 
      2 \[Rho] (a Sin[\[Theta]] Cos[\[Phi]] + 
         b Sin[\[Theta]] Sin[\[Phi]] + c Cos[\[Theta]]) + a^2 + b^2 + 
      c^2 - rn^2 == 
     0 && \[Rho] - \[Lambda] rn Cos[\[Theta] + \[Xi] Cos[\[Phi]]]^2 ==
      0, {\[Rho], \[Theta]}];
(*evaluation*)
In[8]:= ss = soln1[.5, .5, .5, .2, .5, .1, 1]

During evaluation of In[8]:= NSolve::ifun: Inverse functions are being used by NSolve, so some solutions may not be found; use Reduce for complete solution information.

Out[8]= {{\[Rho] -> -2.51366 - 3.61652 I, \[Theta] -> 
   1.8313 - 1.79756 I}, {\[Rho] -> -2.51366 + 3.61652 I, \[Theta] -> 
   1.8313 + 1.79756 I}, {\[Rho] -> -0.146207 - 
    0.0312241 I, \[Theta] -> -1.82025 + 
    0.520643 I}, {\[Rho] -> -0.146207 + 
    0.0312241 I, \[Theta] -> -1.82025 - 0.520643 I}, {\[Rho] -> 
   0.152974, \[Theta] -> 0.785684}, {\[Rho] -> 
   0.311045 - 0.210578 I, \[Theta] -> 
   2.25043 - 0.388346 I}, {\[Rho] -> 
   0.311045 + 0.210578 I, \[Theta] -> 
   2.25043 + 0.388346 I}, {\[Rho] -> 0.417435, \[Theta] -> -0.617471}}
NOTE: The two real solutions are: In[8]:= ss = soln1[.5, .5, .5, .2, .5, .1, 1]

During evaluation of In[8]:= NSolve::ifun: Inverse functions are being used by NSolve, so some solutions may not be found; use Reduce for complete solution information.

Out[8]= {{\[Rho] -> -2.51366 - 3.61652 I, \[Theta] -> 
   1.8313 - 1.79756 I}, {\[Rho] -> -2.51366 + 3.61652 I, \[Theta] -> 
   1.8313 + 1.79756 I}, {\[Rho] -> -0.146207 - 
    0.0312241 I, \[Theta] -> -1.82025 + 
    0.520643 I}, {\[Rho] -> -0.146207 + 
    0.0312241 I, \[Theta] -> -1.82025 - 0.520643 I}, {\[Rho] -> 
   0.152974, \[Theta] -> 0.785684}, {\[Rho] -> 
   0.311045 - 0.210578 I, \[Theta] -> 
   2.25043 - 0.388346 I}, {\[Rho] -> 
   0.311045 + 0.210578 I, \[Theta] -> 
   2.25043 + 0.388346 I}, {\[Rho] -> 0.417435, \[Theta] -> -0.617471}}

NOTE: The two real solutions are:  {\[Rho] -> 
   0.152974, \[Theta] -> 0.785684}
{\[Rho] -> 0.417435, \[Theta] -> -0.617471}
These are the two solution I want.
(2) Using your suggestion (I use soln2 to avoid conflict)
soln2[a_, b_, c_, \[Xi]_, \[Lambda]_, \[Phi]_, rn_] := 
  NSolve[\[Rho]^2 + 
      2 \[Rho] (a Sin[\[Theta]] Cos[\[Phi]] + 
         b Sin[\[Theta]] Sin[\[Phi]] + c Cos[\[Theta]]) + a^2 + b^2 + 
      c^2 - rn^2 == 
     0 && \[Rho] - \[Lambda] rn Cos[\[Theta] + \[Xi] Cos[\[Phi]]]^2 ==
      0, {\[Rho], \[Theta]}, Reals];
EVALUATE]
In[20]:= ss2 = soln2[1/2, 1/2, 1/2, 1/5, 1/2, 1/10, 1]

Out[20]= {{\[Rho] -> 
   ConditionalExpression[
    0.5 Cos[0.199001 + 1. (-0.617471 + 6.28319 C[1])]^2, 
    C[1] \[Element] Integers], \[Theta] -> 
   ConditionalExpression[1. (-0.617471 + 6.28319 C[1]), 
    C[1] \[Element] Integers]}, {\[Rho] -> 
   ConditionalExpression[
    0.5 Cos[0.199001 + 1. (0.785684 + 6.28319 C[1])]^2, 
    C[1] \[Element] Integers], \[Theta] -> 
   ConditionalExpression[1. (0.785684 + 6.28319 C[1]), 
    C[1] \[Element] Integers]}}
Please note: The two real solutions:  {\[Rho] -> 
   0.152974, \[Theta] -> 0.785684}
{\[Rho] -> 0.417435, \[Theta] -> -0.617471}

are no where to be found. Can you illuminate me on this issue? Thank you.

Simply add Real in NSolve:

soln1[a_, b_, c_, \[Xi]_, \[Lambda]_, \[Phi]_, rn_] := 
  NSolve[\[Rho]^2 + 
      2 \[Rho] (a Sin[\[Theta]] Cos[\[Phi]] + 
         b Sin[\[Theta]] Sin[\[Phi]] + c Cos[\[Theta]]) + a^2 + b^2 + 
      c^2 - rn^2 == 
     0 && \[Rho] - \[Lambda] rn Cos[\[Theta] + \[Xi] Cos[\[Phi]]]^2 ==
      0, {\[Rho], \[Theta]}, Reals];
soln1[1/2, 1/2, 1/2, 1/2, 1/2, 1/10, 1]