Home | History | Annotate | Line # | Download | only in src
      1 /* mpc_mul_i -- Multiply a complex number by plus or minus i.
      2 
      3 Copyright (C) 2005, 2009, 2010, 2011, 2012 INRIA
      4 
      5 This file is part of GNU MPC.
      6 
      7 GNU MPC is free software; you can redistribute it and/or modify it under
      8 the terms of the GNU Lesser General Public License as published by the
      9 Free Software Foundation; either version 3 of the License, or (at your
     10 option) any later version.
     11 
     12 GNU MPC is distributed in the hope that it will be useful, but WITHOUT ANY
     13 WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
     14 FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
     15 more details.
     16 
     17 You should have received a copy of the GNU Lesser General Public License
     18 along with this program. If not, see http://www.gnu.org/licenses/ .
     19 */
     20 
     21 #include "mpc-impl.h"
     22 
     23 int
     24 mpc_mul_i (mpc_ptr a, mpc_srcptr b, int sign, mpc_rnd_t rnd)
     25 /* if sign is >= 0, multiply by i, otherwise by -i */
     26 {
     27   int   inex_re, inex_im;
     28   mpfr_t tmp;
     29 
     30   /* Treat the most probable case of compatible precisions first */
     31   if (     MPC_PREC_RE (b) == MPC_PREC_IM (a)
     32         && MPC_PREC_IM (b) == MPC_PREC_RE (a))
     33   {
     34      if (a == b)
     35         mpfr_swap (mpc_realref (a), mpc_imagref (a));
     36      else
     37      {
     38         mpfr_set (mpc_realref (a), mpc_imagref (b), MPFR_RNDN);
     39         mpfr_set (mpc_imagref (a), mpc_realref (b), MPFR_RNDN);
     40      }
     41      if (sign >= 0)
     42         MPFR_CHANGE_SIGN (mpc_realref (a));
     43      else
     44         MPFR_CHANGE_SIGN (mpc_imagref (a));
     45      inex_re = 0;
     46      inex_im = 0;
     47   }
     48   else
     49   {
     50      if (a == b)
     51      {
     52         mpfr_init2 (tmp, MPC_PREC_RE (a));
     53         if (sign >= 0)
     54         {
     55            inex_re = mpfr_neg (tmp, mpc_imagref (b), MPC_RND_RE (rnd));
     56            inex_im = mpfr_set (mpc_imagref (a), mpc_realref (b), MPC_RND_IM (rnd));
     57         }
     58         else
     59         {
     60            inex_re = mpfr_set (tmp, mpc_imagref (b), MPC_RND_RE (rnd));
     61            inex_im = mpfr_neg (mpc_imagref (a), mpc_realref (b), MPC_RND_IM (rnd));
     62         }
     63         mpfr_clear (mpc_realref (a));
     64         mpc_realref (a)[0] = tmp [0];
     65      }
     66      else
     67         if (sign >= 0)
     68         {
     69            inex_re = mpfr_neg (mpc_realref (a), mpc_imagref (b), MPC_RND_RE (rnd));
     70            inex_im = mpfr_set (mpc_imagref (a), mpc_realref (b), MPC_RND_IM (rnd));
     71         }
     72         else
     73         {
     74            inex_re = mpfr_set (mpc_realref (a), mpc_imagref (b), MPC_RND_RE (rnd));
     75            inex_im = mpfr_neg (mpc_imagref (a), mpc_realref (b), MPC_RND_IM (rnd));
     76         }
     77   }
     78 
     79   return MPC_INEX(inex_re, inex_im);
     80 }
     81