/* soft-fp x * y + z as ternary operation.

   Copyright (C) 2006 Free Software Foundation, Inc.

   This file is part of the GNU C Library.

   Contributed by Steven Munroe <sjmunroe@us.ibm.com>, 2006.

   The GNU C Library is free software; you can redistribute it and/or

   modify it under the terms of the GNU Lesser General Public

   License as published by the Free Software Foundation; either

   version 2.1 of the License, or (at your option) any later version.

   The GNU C Library is distributed in the hope that it will be useful,

   but WITHOUT ANY WARRANTY; without even the implied warranty of

   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU

   Lesser General Public License for more details.

   You should have received a copy of the GNU Lesser General Public

   License along with the GNU C Library; if not, write to the Free

   Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA

   02111-1307 USA.  */

#include "soft-fp.h"

#include "double.h"

#include "quad.h"

/* Compute floating point multiply-add with higher (quad) precision.  */

DFtype

__fmadf4 (DFtype a, DFtype b, DFtype c)

{

    FP_DECL_EX;

    FP_DECL_D(A);

    FP_DECL_D(B);

    FP_DECL_D(C);

    FP_DECL_Q(X);

    FP_DECL_Q(Y);

    FP_DECL_Q(Z);

    FP_DECL_Q(U);

    FP_DECL_Q(V);

    FP_DECL_D(R);

    double r;

    long double u, x, y, z;

    FP_INIT_ROUNDMODE;

    FP_UNPACK_RAW_D (A, a);

    FP_UNPACK_RAW_D (B, b);

    FP_UNPACK_RAW_D (C, c);

    /* Extend double to quad.  */

#if (2 * _FP_W_TYPE_SIZE) < _FP_FRACBITS_Q

    FP_EXTEND(Q,D,4,2,X,A);

    FP_EXTEND(Q,D,4,2,Y,B);

    FP_EXTEND(Q,D,4,2,Z,C);

#else

    FP_EXTEND(Q,D,2,1,X,A);

    FP_EXTEND(Q,D,2,1,Y,B);

    FP_EXTEND(Q,D,2,1,Z,C);

#endif

    FP_PACK_RAW_Q(x,X);

    FP_PACK_RAW_Q(y,Y);

    FP_PACK_RAW_Q(z,Z);

    FP_HANDLE_EXCEPTIONS;

    /* Multiply.  */

    FP_INIT_ROUNDMODE;

    FP_UNPACK_Q(X,x);

    FP_UNPACK_Q(Y,y);

    FP_MUL_Q(U,X,Y);

    FP_PACK_Q(u,U);

    FP_HANDLE_EXCEPTIONS;

    /* Subtract.  */

    FP_INIT_ROUNDMODE;

    FP_UNPACK_SEMIRAW_Q(U,u);

    FP_UNPACK_SEMIRAW_Q(Z,z);

    FP_ADD_Q(V,U,Z);

    /* Truncate quad to double.  */

#if (2 * _FP_W_TYPE_SIZE) < _FP_FRACBITS_Q

    V_f[3] &= 0x0007ffff;

    FP_TRUNC(D,Q,2,4,R,V);

#else

    V_f1 &= 0x0007ffffffffffffL;

    FP_TRUNC(D,Q,1,2,R,V);

#endif

    FP_PACK_SEMIRAW_D(r,R);

    FP_HANDLE_EXCEPTIONS;

    return r;

}

/* soft-fp x * y + z as ternary operation.

   Copyright (C) 2006 Free Software Foundation, Inc.

   This file is part of the GNU C Library.

   Contributed by Steven Munroe <sjmunroe@us.ibm.com>, 2006.

   The GNU C Library is free software; you can redistribute it and/or

   modify it under the terms of the GNU Lesser General Public

   License as published by the Free Software Foundation; either

   version 2.1 of the License, or (at your option) any later version.

   The GNU C Library is distributed in the hope that it will be useful,

   but WITHOUT ANY WARRANTY; without even the implied warranty of

   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU

   Lesser General Public License for more details.

   You should have received a copy of the GNU Lesser General Public

   License along with the GNU C Library; if not, write to the Free

   Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA

   02111-1307 USA.  */

#include "soft-fp.h"

#include "single.h"

#include "double.h"

/* Compute floating point multiply-add with higher (double) precision.  */

SFtype

__fmasf4 (SFtype a, SFtype b, SFtype c)

{

    FP_DECL_EX;

    FP_DECL_S(A);

    FP_DECL_S(B);

    FP_DECL_S(C);

    FP_DECL_D(X);

    FP_DECL_D(Y);

    FP_DECL_D(Z);

    FP_DECL_D(U);

    FP_DECL_D(V);

    FP_DECL_S(R);

    float r;

    double u, x, y, z;

    FP_INIT_ROUNDMODE;

    FP_UNPACK_RAW_S (A, a);

    FP_UNPACK_RAW_S (B, b);

    FP_UNPACK_RAW_S (C, c);

    /* Extend single to double.  */

#if _FP_W_TYPE_SIZE < _FP_FRACBITS_D

    FP_EXTEND(D,S,2,1,X,A);

    FP_EXTEND(D,S,2,1,Y,B);

    FP_EXTEND(D,S,2,1,Z,C);

#else

    FP_EXTEND(D,S,1,1,X,A);

    FP_EXTEND(D,S,1,1,Y,B);

    FP_EXTEND(D,S,1,1,Z,C);

#endif

    FP_PACK_RAW_D(x,X);

    FP_PACK_RAW_D(y,Y);

    FP_PACK_RAW_D(z,Z);

    FP_HANDLE_EXCEPTIONS;

    /* Multiply.  */

    FP_INIT_ROUNDMODE;

    FP_UNPACK_D(X,x);

    FP_UNPACK_D(Y,y);

    FP_MUL_D(U,X,Y);

    FP_PACK_D(u,U);

    FP_HANDLE_EXCEPTIONS;

    /* Subtract.  */

    FP_INIT_ROUNDMODE;

    FP_UNPACK_SEMIRAW_D(U,u);

    FP_UNPACK_SEMIRAW_D(Z,z);

    FP_ADD_D(V,U,Z);

    /* Truncate double to single.  */

#if FP_W_TYPE_SIZE < _FP_FRACBITS_D

    V_f1 &= 0x007fffffL;

    FP_TRUNC(S,D,1,2,R,V);

#else

    V_f &= 0x007fffffffffffffL;

    FP_TRUNC(S,D,1,1,R,V);

#endif

    FP_PACK_SEMIRAW_S(r,R);

    FP_HANDLE_EXCEPTIONS;

    return r;

}