src/sqrt.c: Returns value according to ISO C99 specifications for the csqrt function.

tests/tsqrt.c: Add tests for special values.



git-svn-id: svn://scm.gforge.inria.fr/svn/mpc/trunk@156 211d60ee-9f03-0410-a15a-8952a2c7a4e4

1 files changed, 649 insertions, 1 deletions

diff --git a/tests/tsqrt.c b/tests/tsqrt.c
index 9e1fbbd..9bb0067 100644
--- a/tests/tsqrt.c
+++ b/tests/tsqrt.c
@@ -1,6 +1,6 @@
 /* tsqrt -- test file for mpc_sqrt.
 
-Copyright (C) 2008 INRIA.
+Copyright (C) 2008 Philippe Th\'eveny.
 
 This file is part of the MPC Library.
 
@@ -30,11 +30,659 @@ MA 02111-1307, USA. */
 #define TEST_FUNCTION mpc_sqrt
 #include "tgeneric.c"
 
+static void
+test_failed (mpc_t op, mpc_t get, mpc_t expected)
+{
+  printf ("mpc_sqrt(op) failed\n with ");
+  OUT (op);
+  printf ("     ");
+  OUT (get);
+  OUT (expected);
+  exit (1);
+}
+
+/* check special values as defined in C99 standard */
+static void
+special ()
+{
+  mpc_t z;
+  mpc_t t;
+
+  mpc_init (z);
+  mpc_init (t);
+
+  /* sqrt(+Inf +i*Inf) = +Inf +i*Inf */
+  mpfr_set_inf (MPC_RE (z), +1);
+  mpfr_set_inf (MPC_IM (z), +1);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_inf_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_inf_p (MPC_IM (t)) || mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_inf (MPC_RE(c99), +1);
+      mpfr_set_inf (MPC_IM(c99), +1);
+      test_failed (z, t, c99);
+    }
+
+  /* sqrt(+Inf -i*Inf) =  +Inf -i*Inf */
+  mpc_conj (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_inf_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_inf_p (MPC_IM (t)) || !mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_inf (MPC_RE(c99), +1);
+      mpfr_set_inf (MPC_IM(c99), -1);
+      test_failed (z, t, c99);
+    }
+
+  /* sqrt(-Inf +i*Inf) = +Inf +i*Inf */
+  mpc_neg (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_inf_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_inf_p (MPC_IM (t)) || mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_inf (MPC_RE(c99), +1);
+      mpfr_set_inf (MPC_IM(c99), +1);
+      test_failed (z, t, c99);
+    }
+
+  /* sqrt(-Inf -i*Inf) = +Inf -i*Inf */
+  mpc_conj (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_inf_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_inf_p (MPC_IM (t)) || !mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_inf (MPC_RE(c99), +1);
+      mpfr_set_inf (MPC_IM(c99), -1);
+      test_failed (z, t, c99);
+    }
+
+  /* tan (NaN -i*Inf) = +Inf -i*Inf */
+  mpfr_set_nan (MPC_RE (z));
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_inf_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_inf_p (MPC_IM (t)) || !mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_inf (MPC_RE(c99), +1);
+      mpfr_set_inf (MPC_IM(c99), -1);
+      test_failed (z, t, c99);
+    }
+
+  /* tan (NaN +i*Inf) = +Inf +i*Inf */
+  mpc_conj (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_inf_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_inf_p (MPC_IM (t)) || mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_inf (MPC_RE(c99), +1);
+      mpfr_set_inf (MPC_IM(c99), +1);
+      test_failed (z, t, c99);
+    }
+
+  /* sqrt(1 +i*Inf) = +Inf +i*Inf */
+  mpfr_set_ui (MPC_RE (z), 1, GMP_RNDN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_inf_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_inf_p (MPC_IM (t)) || mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_inf (MPC_RE(c99), +1);
+      mpfr_set_inf (MPC_IM(c99), +1);
+      test_failed (z, t, c99);
+    }
+
+  /* sqrt(1 -i*Inf) = +Inf -i*Inf */
+  mpc_conj (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_inf_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_inf_p (MPC_IM (t)) || !mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_inf (MPC_RE(c99), +1);
+      mpfr_set_inf (MPC_IM(c99), -1);
+      test_failed (z, t, c99);
+    }
+
+  /* sqrt(-1 +i*Inf) = +Inf +i*Inf */
+  mpc_neg (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_inf_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_inf_p (MPC_IM (t)) || mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_inf (MPC_RE(c99), +1);
+      mpfr_set_inf (MPC_IM(c99), +1);
+      test_failed (z, t, c99);
+    }
+
+  /* sqrt(-1 -i*Inf) = +Inf -i*Inf */
+  mpc_conj (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_inf_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_inf_p (MPC_IM (t)) || !mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_inf (MPC_RE(c99), +1);
+      mpfr_set_inf (MPC_IM(c99), -1);
+      test_failed (z, t, c99);
+    }
+
+  /*  sqrt(+0 -i*Inf) = +Inf -i*Inf */
+  mpfr_set_ui (MPC_RE (z), 0, GMP_RNDN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_inf_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_inf_p (MPC_IM (t)) || !mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_inf (MPC_RE(c99), +1);
+      mpfr_set_inf (MPC_IM(c99), -1);
+      test_failed (z, t, c99);
+    }
+
+  /*  sqrt(+0 +i*Inf) = +Inf +i*Inf */
+  mpc_conj (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_inf_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_inf_p (MPC_IM (t)) || mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_inf (MPC_RE(c99), +1);
+      mpfr_set_inf (MPC_IM(c99), +1);
+      test_failed (z, t, c99);
+    }
+
+  /*  sqrt(-0 -i*Inf) = +Inf -i*Inf */
+  mpc_neg (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_inf_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_inf_p (MPC_IM (t)) || !mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_inf (MPC_RE(c99), +1);
+      mpfr_set_inf (MPC_IM(c99), -1);
+      test_failed (z, t, c99);
+    }
+
+  /*  sqrt(-0 +i*Inf) = +Inf +i*Inf */
+  mpc_conj (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_inf_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_inf_p (MPC_IM (t)) || mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_inf (MPC_RE(c99), +1);
+      mpfr_set_inf (MPC_IM(c99), +1);
+      test_failed (z, t, c99);
+    }
+
+  /* sqrt(-0 +0*i) = +0 +0*i */
+  mpfr_set_ui (MPC_IM (z), 0, GMP_RNDN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_zero_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_zero_p (MPC_IM (t)) || mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpc_set_ui_ui (c99, 0, 0, MPC_RNDNN);
+      test_failed (z, t, c99);
+    }
+
+  /* sqrt(-0 -0*i) = +0 -0*i */
+  mpc_conj (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_zero_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_zero_p (MPC_IM (t)) || !mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpc_set_ui_ui (c99, 0, 0, MPC_RNDNN);
+      mpc_conj (c99, c99, MPC_RNDNN);
+      test_failed (z, t, c99);
+    }
+
+  /* sqrt(+0 +0*i) = +0 +0*i */
+  mpc_neg (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_zero_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_zero_p (MPC_IM (t)) || mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpc_set_ui_ui (c99, 0, 0, MPC_RNDNN);
+      test_failed (z, t, c99);
+    }
+
+  /* sqrt(+0 -0*i) = +0 -0*i */
+  mpc_conj (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_zero_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_zero_p (MPC_IM (t)) || !mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpc_set_ui_ui (c99, 0, 0, MPC_RNDNN);
+      mpc_conj (c99, c99, MPC_RNDNN);
+      test_failed (z, t, c99);
+    }
+
+  /* sqrt(x +i*NaN) = NaN+i*NaN, when x is finite */
+  mpfr_set_ui (MPC_RE (z), 1, GMP_RNDN);
+  mpfr_set_nan (MPC_IM (z));
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_nan_p (MPC_RE (t)) || !mpfr_nan_p (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_nan (MPC_RE (c99));
+      mpfr_set_nan (MPC_IM (c99));
+      test_failed (z, t, c99);
+    }
+
+  mpc_neg (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_nan_p (MPC_RE (t)) || !mpfr_nan_p (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_nan (MPC_RE (c99));
+      mpfr_set_nan (MPC_IM (c99));
+      test_failed (z, t, c99);
+    }
+
+  mpfr_set_ui (MPC_RE (z), 0, GMP_RNDN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_nan_p (MPC_RE (t)) || !mpfr_nan_p (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_nan (MPC_RE (c99));
+      mpfr_set_nan (MPC_IM (c99));
+      test_failed (z, t, c99);
+    }
+
+  mpc_neg (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_nan_p (MPC_RE (t)) || !mpfr_nan_p (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_nan (MPC_RE (c99));
+      mpfr_set_nan (MPC_IM (c99));
+      test_failed (z, t, c99);
+    }
+
+  /* sqrt(NaN +i*y) = NaN+i*NaN, when y is finite */
+  mpfr_set_nan (MPC_RE (z));
+  mpfr_set_ui (MPC_IM (z), 1, GMP_RNDN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_nan_p (MPC_RE (t)) || !mpfr_nan_p (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_nan (MPC_RE (c99));
+      mpfr_set_nan (MPC_IM (c99));
+      test_failed (z, t, c99);
+    }
+
+  mpc_conj (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_nan_p (MPC_RE (t)) || !mpfr_nan_p (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_nan (MPC_RE (c99));
+      mpfr_set_nan (MPC_IM (c99));
+      test_failed (z, t, c99);
+    }
+
+  mpfr_set_ui (MPC_IM (z), 0, GMP_RNDN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_nan_p (MPC_RE (t)) || !mpfr_nan_p (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_nan (MPC_RE (c99));
+      mpfr_set_nan (MPC_IM (c99));
+      test_failed (z, t, c99);
+    }
+
+  mpc_conj (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_nan_p (MPC_RE (t)) || !mpfr_nan_p (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_nan (MPC_RE (c99));
+      mpfr_set_nan (MPC_IM (c99));
+      test_failed (z, t, c99);
+    }
+
+  /* sqrt(Nan +i*NaN) = NaN +i*NaN */
+  mpfr_set_nan (MPC_IM (z));
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_nan_p (MPC_RE (t)) || !mpfr_nan_p (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_nan (MPC_RE (c99));
+      mpfr_set_nan (MPC_IM (c99));
+      test_failed (z, t, c99);
+    }
+
+  /* sqrt(-Inf +i*NaN) = NaN +/-i*Inf */
+  mpfr_set_inf (MPC_RE (z), -1);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_nan_p (MPC_RE (t)) || !mpfr_inf_p (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_nan (MPC_RE (c99));
+      mpfr_set_inf (MPC_IM (c99), -1);
+      test_failed (z, t, c99);
+    }
+
+  /* sqrt(+Inf +i*NaN) = +Inf +i*NaN */
+  mpfr_set_inf (MPC_RE (z), +1);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_inf_p (MPC_RE (t)) || !mpfr_nan_p (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_inf (MPC_RE (c99), +1);
+      mpfr_set_nan (MPC_IM (c99));
+      test_failed (z, t, c99);
+    }
+
+  /* sqrt(+Inf +i*y) = +Inf +i*0, when y is a finite non negative number. */
+  /* sqrt(+Inf +i*y) = +Inf -i*0, when y is a finite non positive number. */
+  mpfr_set_ui (MPC_IM (z), 1, GMP_RNDN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_inf_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_zero_p (MPC_IM (t)) || mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_inf (MPC_RE (c99), +1);
+      mpfr_set_ui (MPC_IM (c99), 0, GMP_RNDN);
+      test_failed (z, t, c99);
+    }
+
+  mpc_conj (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_inf_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_zero_p (MPC_IM (t)) || !mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_inf (MPC_RE (c99), +1);
+      mpfr_set_ui (MPC_IM (c99), 0, GMP_RNDN);
+      mpc_conj (c99, c99, MPC_RNDNN);
+      test_failed (z, t, c99);
+    }
+
+  mpfr_set_ui (MPC_IM (z), 0, GMP_RNDN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_inf_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_zero_p (MPC_IM (t)) || mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_inf (MPC_RE (c99), +1);
+      mpfr_set_ui (MPC_IM (c99), 0, GMP_RNDN);
+      test_failed (z, t, c99);
+    }
+
+  mpc_conj (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_inf_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_zero_p (MPC_IM (t)) || !mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_inf (MPC_RE (c99), +1);
+      mpfr_set_ui (MPC_IM (c99), 0, GMP_RNDN);
+      mpc_conj (c99, c99, MPC_RNDNN);
+      test_failed (z, t, c99);
+    }
+
+  /* sqrt(-Inf +i*y) = +0 -i*Inf, when y is a finite non positive number. */
+  /* sqrt(-Inf +i*y) = +0 +i*Inf, when y is a finite non negative number. */
+  mpc_neg (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_zero_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_inf_p (MPC_IM (t)) || mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_ui (MPC_RE (c99), 0, GMP_RNDN);
+      mpfr_set_inf (MPC_IM (c99), +1);
+      test_failed (z, t, c99);
+    }
+
+  mpc_conj (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_zero_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_inf_p (MPC_IM (t)) || !mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_ui (MPC_RE (c99), 0, GMP_RNDN);
+      mpfr_set_inf (MPC_IM (c99), -1);
+      test_failed (z, t, c99);
+    }
+
+  mpfr_set_ui (MPC_IM (z), 1, GMP_RNDN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_zero_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_inf_p (MPC_IM (t)) || mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_ui (MPC_RE (c99), 0, GMP_RNDN);
+      mpfr_set_inf (MPC_IM (c99), +1);
+      test_failed (z, t, c99);
+    }
+
+  mpc_conj (z, z, MPC_RNDNN);
+  mpc_sqrt (t, z, MPC_RNDNN);
+  if (!mpfr_zero_p (MPC_RE (t)) || mpfr_signbit (MPC_RE (t))
+      || !mpfr_inf_p (MPC_IM (t)) || !mpfr_signbit (MPC_IM (t)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_ui (MPC_RE (c99), 0, GMP_RNDN);
+      mpfr_set_inf (MPC_IM (c99), -1);
+      test_failed (z, t, c99);
+    }
+
+  mpc_clear (t);
+  mpc_clear (z);
+}
+
+static void
+pure_real_argument ()
+{
+  /* sqrt(x +i*0) = sqrt(x) +i*0, when x>0 */
+  /* sqrt(x -i*0) = sqrt(x) -i*0, when x>0 */
+  /* sqrt(x +i*0) = +0 +i*sqrt(-x) +i*0, when x<0 */
+  /* sqrt(x -i*0) = +0 -i*sqrt(-x) +i*0, when x<0 */
+  mpfr_t x;
+  mpfr_t sqrt_x;
+  mpc_t z;
+  mpc_t sqrt_z;
+
+  mpfr_init (x);
+  mpfr_init (sqrt_x);
+  mpc_init (z);
+  mpc_init (sqrt_z);
+
+  /* sqrt(2 +i*0) = sqrt(2) +i*0 */
+  mpc_set_ui_ui (z, 2, 0, MPC_RNDNN);
+  mpfr_set_ui (x, 2, GMP_RNDN);
+  mpfr_sqrt (sqrt_x, x, GMP_RNDN);
+  mpc_sqrt (sqrt_z, z, MPC_RNDNN);
+  if (mpfr_cmp (MPC_RE (sqrt_z), sqrt_x) != 0
+      || !mpfr_zero_p (MPC_IM (sqrt_z)) || mpfr_signbit (MPC_IM (sqrt_z)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set (MPC_RE (c99), sqrt_x, GMP_RNDN);
+      mpfr_set_ui (MPC_IM (c99), 0, GMP_RNDN);
+      test_failed (z, sqrt_z, c99);
+    }
+
+  /* sqrt(2 -i*0) = sqrt(2) -i*0 */
+  mpc_conj (z, z, MPC_RNDNN);
+  mpc_sqrt (sqrt_z, z, MPC_RNDNN);
+  if (mpfr_cmp (MPC_RE (sqrt_z), sqrt_x) != 0
+      || !mpfr_zero_p (MPC_IM (sqrt_z)) || !mpfr_signbit (MPC_IM (sqrt_z)))
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set (MPC_RE (c99), sqrt_x, GMP_RNDN);
+      mpfr_set_ui (MPC_IM (c99), 0, GMP_RNDN);
+      mpc_conj (c99, c99, MPC_RNDNN);
+      test_failed (z, sqrt_z, c99);
+    }
+
+  /* sqrt(-2 +i*0) = +0 +i*sqrt(2) */
+  mpc_neg (z, z, MPC_RNDNN);
+  mpc_sqrt (sqrt_z, z, MPC_RNDNN);
+  if (!mpfr_zero_p (MPC_RE (sqrt_z)) || mpfr_signbit (MPC_RE (sqrt_z))
+      || mpfr_cmp (MPC_IM (sqrt_z), sqrt_x) != 0)
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_ui (MPC_RE (c99), 0, GMP_RNDN);
+      mpfr_set (MPC_IM (c99), sqrt_x, GMP_RNDN);
+      test_failed (z, sqrt_z, c99);
+    }
+
+  /* sqrt(-2 -i*0) = +0 -i*sqrt(2) */
+  mpfr_neg (sqrt_x, sqrt_x, GMP_RNDN);
+  mpc_conj (z, z, MPC_RNDNN);
+  mpc_sqrt (sqrt_z, z, MPC_RNDNN);
+  if (!mpfr_zero_p (MPC_RE (sqrt_z)) || mpfr_signbit (MPC_RE (sqrt_z))
+      || mpfr_cmp (MPC_IM (sqrt_z), sqrt_x) != 0)
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set_ui (MPC_RE (c99), 0, GMP_RNDN);
+      mpfr_set (MPC_IM (c99), sqrt_x, GMP_RNDN);
+      test_failed (z, sqrt_z, c99);
+    }
+
+  mpc_clear (sqrt_z);
+  mpc_clear (z);
+  mpfr_clear (sqrt_x);
+  mpfr_clear (x);
+}
+
+static void
+pure_imaginary_argument ()
+{
+  /* sqrt(+/-0 +i*y) = sqrt(y/2) * (1 +i), when y >0     */
+  /* sqrt(+/-0 +i*y) = sqrt(-y/2) * (-1 + i), when y < 0 */
+  /* here, b = sqrt(y/2) */
+  mpfr_t b;
+  mpfr_t minus_b;
+  mpc_t z;
+  mpc_t sqrt_z;
+
+  mpfr_init (b);
+  mpfr_init (minus_b);
+  mpc_init (z);
+  mpc_init (sqrt_z);
+
+  /* sqrt(+0 +i*4) = sqrt(2) + i* sqrt(2) */
+  mpfr_set_ui (b, 2, GMP_RNDN);
+  mpfr_sqrt (b, b, GMP_RNDN);
+  mpc_set_ui_ui (z, 0, 4, MPC_RNDNN);
+  mpc_sqrt (sqrt_z, z, MPC_RNDNN);
+  if (mpfr_cmp (MPC_RE (sqrt_z), b) != 0
+      || mpfr_cmp (MPC_IM (sqrt_z), b) != 0)
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set (MPC_RE (c99), b, GMP_RNDN);
+      mpfr_set (MPC_IM (c99), b, GMP_RNDN);
+      test_failed (z, sqrt_z, c99);
+    }
+
+  /* sqrt(+0 -i*4) = -sqrt(2) + i* sqrt(2) */
+  mpfr_neg (minus_b, b, GMP_RNDN);
+  mpc_conj (z, z, MPC_RNDNN);
+  mpc_sqrt (sqrt_z, z, MPC_RNDNN);
+  if (mpfr_cmp (MPC_RE (sqrt_z), b) != 0
+      || mpfr_cmp (MPC_IM (sqrt_z), minus_b) != 0)
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set (MPC_RE (c99), b, GMP_RNDN);
+      mpfr_set (MPC_IM (c99), minus_b, GMP_RNDN);
+      test_failed (z, sqrt_z, c99);
+    }
+
+  /* sqrt(-0 +i*4) = sqrt(2) + i* sqrt(2) */
+  mpc_neg (z, z, MPC_RNDNN);
+  mpc_sqrt (sqrt_z, z, MPC_RNDNN);
+  if (mpfr_cmp (MPC_RE (sqrt_z), b) != 0
+      || mpfr_cmp (MPC_IM (sqrt_z), b) != 0)
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set (MPC_RE (c99), b, GMP_RNDN);
+      mpfr_set (MPC_IM (c99), b, GMP_RNDN);
+      test_failed (z, sqrt_z, c99);
+    }
+
+  /* sqrt(-0 -i*4) = -sqrt(2) + i* sqrt(2) */
+  mpc_conj (z, z, MPC_RNDNN);
+  mpc_sqrt (sqrt_z, z, MPC_RNDNN);
+  if (mpfr_cmp (MPC_RE (sqrt_z), b) != 0
+      || mpfr_cmp (MPC_IM (sqrt_z), minus_b) != 0)
+    {
+      mpc_t c99;
+      mpc_init (c99);
+      mpfr_set (MPC_RE (c99), b, GMP_RNDN);
+      mpfr_set (MPC_IM (c99), minus_b, GMP_RNDN);
+      test_failed (z, sqrt_z, c99);
+    }
+
+  mpc_clear (sqrt_z);
+  mpc_clear (z);
+  mpfr_clear (minus_b);
+  mpfr_clear (b);
+}
+
 int
 main()
 {
   test_start ();
 
+  special ();
+  pure_real_argument ();
+  pure_imaginary_argument ();
+
   tgeneric ();
 
   test_end ();