static char rcsid[] = "$Id: F-opt-expv.c,v 1.21 2002/01/13 15:52:26 a-hasega Exp $";
/* 
 * $RWC_Release: Omni-1.6 $
 * $RWC_Copyright:
 *  Omni Compiler Software Version 1.5-1.6
 *  Copyright (C) 2002 PC Cluster Consortium
 *  
 *  This software is free software; you can redistribute it and/or modify
 *  it under the terms of the GNU Lesser General Public License version
 *  2.1 published by the Free Software Foundation.
 *  
 *  Omni Compiler Software Version 1.0-1.4
 *  Copyright (C) 1999, 2000, 2001.
 *   Tsukuba Research Center, Real World Computing Partnership, Japan.
 *  
 *  Please check the Copyright and License information in the files named
 *  COPYRIGHT and LICENSE under the top  directory of the Omni Compiler
 *  Software release kit.
 *  
 *  
 *  $
 */
#include "F-front.h"
#include <math.h>
#define REDUCE_FLOAT 1

#define EXPV_IS_ZERO(v)	\
  (EXPV_CODE(v) == INT_CONSTANT && EXPV_INT_VALUE(v) == 0)
#define EXPV_IS_ONE(v)	\
  (EXPV_CODE(v) == INT_CONSTANT && EXPV_INT_VALUE(v) == 1)
#define EXPV_IS_MINUS_ONE(v)	\
  (EXPV_CODE(v) == INT_CONSTANT && EXPV_INT_VALUE(v) == -1)
#define EXPV_IS_ZERO_FLOAT(v)	\
  (EXPV_CODE(v) == FLOAT_CONSTANT && EXPV_FLOAT_VALUE(v) == 0.0)
#define EXPV_IS_ONE_FLOAT(v)	\
  (EXPV_CODE(v) == FLOAT_CONSTANT && EXPV_FLOAT_VALUE(v) == 1.0)


expv
expv_numeric_const_reduce(left, right, code, v)
     expv left;
     expv right;
     enum expr_code code;
     expv v;
{
    BASIC_DATA_TYPE bTyp;
    TYPE_DESC tp = NULL;

    expv nL = NULL;
    expv nR = NULL;

    if (v != NULL) {
	tp = EXPV_TYPE(v);
    }
    if (tp != NULL) {
	bTyp = TYPE_BASIC_TYPE(tp);
    } else {
	if (right != NULL) {
	    tp = max_type(EXPV_TYPE(left), EXPV_TYPE(right));
	} else {
	    tp = EXPV_TYPE(left);
	}
	bTyp = TYPE_BASIC_TYPE(tp);
    }

    if (IS_NUMERIC_CONST_V(left)) {
	nL = expv_reduce_conv_const(tp, left);
    }
    if (right != NULL) {
	if (IS_NUMERIC_CONST_V(right)) {
	    nR = expv_reduce_conv_const(tp, right);
	}
    }

    switch (bTyp) {
	case TYPE_SHORT:
	case TYPE_INT: {
	    int i;
	    
	    switch (code) {
		case MUL_EXPR: {
		    i = EXPV_INT_VALUE(nL) * EXPV_INT_VALUE(nR);
		    break;
		}
		case DIV_EXPR: {
		    if (EXPV_INT_VALUE(nR) == 0) {
			error_at_node(v, "divide by zero");
			goto NonReducedReturn;
		    }
		    i = EXPV_INT_VALUE(nL) / EXPV_INT_VALUE(nR);
		    break;
		}
		case PLUS_EXPR: {
		    i = EXPV_INT_VALUE(nL) + EXPV_INT_VALUE(nR);
		    break;
		}
		case MINUS_EXPR: {
		    i = EXPV_INT_VALUE(nL) - EXPV_INT_VALUE(nR);
		    break;
		}
		case UNARY_MINUS_EXPR: {
		    i = -EXPV_INT_VALUE(nL);
		    break;
		}
		case MOD_EXPR: {
		    if (EXPV_INT_VALUE(nR) == 0) {
			error_at_node((expr)v, "mod by zero");
			goto NonReducedReturn;
		    }
		    i = EXPV_INT_VALUE(nL) % EXPV_INT_VALUE(nR);
		    break;
		}
		case LSHIFT_EXPR: {
		    if (EXPV_INT_VALUE(nR) == 0) {
			warning_at_node((expr)v, "left shift by zero");
		    } else if (EXPV_INT_VALUE(nR) >= 32) {
			warning_at_node((expr)v, "left shift by >= 32");
		    }
		    i = (EXPV_INT_VALUE(nL) << EXPV_INT_VALUE(nR));
		    break;
		}
		case RSHIFT_EXPR: {
		    if (EXPV_INT_VALUE(nR) == 0) {
			warning_at_node((expr)v, "right shift by zero");
		    } else if (EXPV_INT_VALUE(nR) >= 32) {
			warning_at_node((expr)v, "right shift by >= 32");
		    }
		    i = (EXPV_INT_VALUE(nL) >> EXPV_INT_VALUE(nR));
		    break;
		}
		case BIT_AND_EXPR: {
		    i = (EXPV_INT_VALUE(nL) & EXPV_INT_VALUE(nR));
		    break;
		}
		case BIT_OR_EXPR: {
		    i = (EXPV_INT_VALUE(nL) | EXPV_INT_VALUE(nR));
		    break;
		}
		case BIT_XOR_EXPR: {
		    i = (EXPV_INT_VALUE(nL) ^ EXPV_INT_VALUE(nR));
		    break;
		}
		case BIT_NOT_EXPR: {
		    i = ~(EXPV_INT_VALUE(nL));
		    break;
		}
		case LOG_EQ_EXPR: {
		    i = (EXPV_INT_VALUE(nL) == EXPV_INT_VALUE(nR));
		    break;
		}
		case LOG_NEQ_EXPR: {
		    i = (EXPV_INT_VALUE(nL) != EXPV_INT_VALUE(nR));
		    break;
		}
		case LOG_GE_EXPR: {
		    i = (EXPV_INT_VALUE(nL) >= EXPV_INT_VALUE(nR));
		    break;
		}
		case LOG_GT_EXPR: {
		    i = (EXPV_INT_VALUE(nL) > EXPV_INT_VALUE(nR));
		    break;
		}
		case LOG_LE_EXPR: {
		    i = (EXPV_INT_VALUE(nL) <= EXPV_INT_VALUE(nR));
		    break;
		}
		case LOG_LT_EXPR: {
		    i = (EXPV_INT_VALUE(nL) < EXPV_INT_VALUE(nR));
		    break;
		}
		case LOG_AND_EXPR: {
		    i = (EXPV_INT_VALUE(nL) && EXPV_INT_VALUE(nR));
		    break;
		}
		case LOG_OR_EXPR: {
		    i = (EXPV_INT_VALUE(nL) || EXPV_INT_VALUE(nR));
		    break;
		}
		case LOG_NOT_EXPR: {
		    i = (!(EXPV_INT_VALUE(nL)));
		    break;
		}
		default: {
		    goto NonReducedReturn;
		}
	    }
	    return expv_int_term(INT_CONSTANT, tp, i);
	}

#ifdef HAS_INT64
	case TYPE_LONGLONG: {
	    _omInt64_t i64 = (_omInt64_t)0;
	    _omInt64_t iL64 = EXPV_INT64_VALUE(nL);
	    _omInt64_t iR64 = (nR != NULL) ? (EXPV_INT64_VALUE(nR)) : 0;
	    expv ret = expv_longlong_term(LONGLONG_CONSTANT, tp, 0, 0);

	    switch (code) {
		case MUL_EXPR: {
		    i64 = iL64 * iR64;
		    break;
		}
		case DIV_EXPR: {
		    if (iR64 == 0) {
			error_at_node(v, "divide by zero");
			goto NonReducedReturn;
		    }
		    i64 = iL64 / iR64;
		    break;
		}
		case PLUS_EXPR: {
		    i64 = iL64 + iR64;
		    break;
		}
		case MINUS_EXPR: {
		    i64 = iL64 - iR64;
		    break;
		}
		case UNARY_MINUS_EXPR: {
		    i64 = -iL64;
		    break;
		}
		case MOD_EXPR: {
		    if (iR64 == 0) {
			error_at_node((expr)v, "mod by zero");
			goto NonReducedReturn;
		    }
		    i64 = iL64 % iR64;
		    break;
		}
		case LSHIFT_EXPR: {
		    if (iR64 == 0) {
			warning_at_node((expr)v, "left shift by zero");
		    } else if (iR64 >= 64) {
			warning_at_node((expr)v, "left shift by >= 64");
		    }
		    i64 = (iL64 << iR64);
		    break;
		}
		case RSHIFT_EXPR: {
		    if (iR64 == 0) {
			warning_at_node((expr)v, "right shift by zero");
		    } else if (iR64 >= 64) {
			warning_at_node((expr)v, "right shift by >= 64");
		    }
		    iR64 = (iL64 >> iR64);
		    break;
		}
		case BIT_AND_EXPR: {
		    i64 = (iL64 & iR64);
		    break;
		}
		case BIT_OR_EXPR: {
		    i64 = (iL64 | iR64);
		    break;
		}
		case BIT_XOR_EXPR: {
		    i64 = (iL64 ^ iR64);
		    break;
		}
		case BIT_NOT_EXPR: {
		    i64 = ~(iL64);
		    break;
		}
		case LOG_EQ_EXPR: {
		    i64 = (iL64 == iR64);
		    break;
		}
		case LOG_NEQ_EXPR: {
		    i64 = (iL64 != iR64);
		    break;
		}
		case LOG_GE_EXPR: {
		    i64 = (iL64 >= iR64);
		    break;
		}
		case LOG_GT_EXPR: {
		    i64 = (iL64 > iR64);
		    break;
		}
		case LOG_LE_EXPR: {
		    i64 = (iL64 <= iR64);
		    break;
		}
		case LOG_LT_EXPR: {
		    i64 = (iL64 < iR64);
		    break;
		}
		case LOG_AND_EXPR: {
		    i64 = (iL64 && iR64);
		    break;
		}
		case LOG_OR_EXPR: {
		    i64 = (iL64 || iR64);
		    break;
		}
		case LOG_NOT_EXPR: {
		    i64 = (!(iL64));
		    break;
		}
		default: {
		    goto NonReducedReturn;
		}
	    }
	    EXPV_INT64_SAVE_VALUE(ret, i64);
	    return ret;
	}
#endif /* HAS_INT64 */

	default: {
	    goto NonReducedReturn;
	}
    }

    NonReducedReturn:
    if (code == EXPV_CODE(v) && tp == EXPV_TYPE(v) &&
	left == EXPV_LEFT(v) && right == EXPV_RIGHT(v)) {
	return v;
    } else {
	return expv_cons(code, tp, left, right);
    }
}


/* 
 * optimize expression value
 */
expv
expv_reduce(v)
     expv v;
{
    enum expr_code code;
    TYPE_DESC tp;
    expv left,right,rv;
    list lp;
    
    if (v == NULL) return(v);	/* error recovery */

    code = EXPV_CODE(v);
    tp = EXPV_TYPE(v);

    /* check for terminal */
    if (EXPR_CODE_IS_TERMINAL_OR_CONST(code)) {
	return v;
    }
    
    if(EXPR_CODE_IS_LIST(code)){
	/* expand list. */
	FOR_ITEMS_IN_LIST(lp,v)
	    LIST_ITEM(lp) = expv_reduce(LIST_ITEM(lp));
	return v;
    }

    /* internal node */
    left = expv_reduce(EXPV_LEFT(v));
    right = EXPV_RIGHT(v);
    if(right != NULL) right = expv_reduce(right);

    /* constant folding */
    switch(code){
    case CAST_EXPR:
	if (IS_NUMERIC_CONST_V(left)) {
	    if ((rv = expv_reduce_conv_const(tp, left)) != NULL) {
		return rv;
	    }
	}
	break;

    case MUL_EXPR:
	if(EXPV_IS_ZERO(left) || EXPV_IS_ZERO(right))
	  return(expv_constant_0); /* x*0 = 0 */
	if(EXPV_IS_ZERO_FLOAT(left) || EXPV_IS_ZERO_FLOAT(right))
	  return(expv_float_0);	/* x*0 = 0 */
	if(EXPV_IS_ONE(left) || EXPV_IS_ONE_FLOAT(left)) 
	  return(right);	/* x*1 = x */
	if(EXPV_IS_ONE(right) || EXPV_IS_ONE_FLOAT(right)) 
	  return(left);		/* 1*x = x */

	if(EXPV_CODE(left) == INT_CONSTANT &&
	   EXPV_CODE(right) == INT_CONSTANT){
	    return(expv_int_term(INT_CONSTANT,EXPV_TYPE(v),
				 EXPV_INT_VALUE(left)*EXPV_INT_VALUE(right)));
	}
#ifdef REDUCE_FLOAT
	if(EXPV_CODE(left) == FLOAT_CONSTANT &&
	   EXPV_CODE(right) == FLOAT_CONSTANT)
	  return(expv_float_term(FLOAT_CONSTANT,EXPV_TYPE(v),
				 EXPV_FLOAT_VALUE(left)*EXPV_FLOAT_VALUE(right)));
#endif
	if (IS_NUMERIC_CONST_V(left) &&
	    IS_NUMERIC_CONST_V(right)) {
	    return expv_numeric_const_reduce(left, right, code, v);
	}
	break;

    case DIV_EXPR:
	if(EXPV_IS_ZERO(left)) return(expv_constant_0); /* 0/x = 0 */
	if(EXPV_IS_ZERO_FLOAT(left)) return(expv_float_0); /* 0/x = 0 */
	if(EXPV_IS_ZERO(right) || EXPV_IS_ZERO_FLOAT(left)) /* x/0 = error */
	  {
	      error_at_node((expr)v, "divide by zero");
	      return(v);
	  }
	if(EXPV_IS_ONE(right) || EXPV_IS_ONE_FLOAT(right))
	  return(left);		/* x/1 = x */

	if(EXPV_CODE(left) == INT_CONSTANT && 
	   EXPV_CODE(right) == INT_CONSTANT)
	    return(expv_int_term(INT_CONSTANT,EXPV_TYPE(v),
				 EXPV_INT_VALUE(left)/EXPV_INT_VALUE(right)));
#ifdef REDUCE_FLOAT
	if(EXPV_CODE(left) == FLOAT_CONSTANT &&
	   EXPV_CODE(right) == FLOAT_CONSTANT)
	  return(expv_float_term(FLOAT_CONSTANT,EXPV_TYPE(v),
				 EXPV_FLOAT_VALUE(left)/EXPV_FLOAT_VALUE(right)));
#endif
	if (IS_NUMERIC_CONST_V(left) &&
	    IS_NUMERIC_CONST_V(right)) {
	    return expv_numeric_const_reduce(left, right, code, v);
	}
	break;

    case PLUS_EXPR:
	if(EXPV_IS_ZERO(left) || EXPV_IS_ZERO_FLOAT(left)) 
	  return(right);	/* 0 + x = x */
	if(EXPV_IS_ZERO(right) || EXPV_IS_ZERO_FLOAT(right)) 
	  return(left);		/* x + 0 = x */

	if(EXPV_CODE(left) == INT_CONSTANT && 
	   EXPV_CODE(right) == INT_CONSTANT)
	  return(expv_int_term(INT_CONSTANT,EXPV_TYPE(v),
			       EXPV_INT_VALUE(left)+EXPV_INT_VALUE(right)));
#ifdef REDUCE_FLOAT
	if(EXPV_CODE(left) == FLOAT_CONSTANT && 
	   EXPV_CODE(right) == FLOAT_CONSTANT)
	  return(expv_float_term(FLOAT_CONSTANT,EXPV_TYPE(v),
				 EXPV_FLOAT_VALUE(left)+EXPV_FLOAT_VALUE(right)));
#endif
	/*  (PLUS (PLUS x c1) c2) => (PLUS x c1+c2) */
	if(EXPV_CODE(right) == INT_CONSTANT  &&
	   EXPV_CODE(left) == PLUS_EXPR &&
	   EXPV_CODE((expv)EXPV_RIGHT(left)) == INT_CONSTANT)
	  return(expv_cons(code,tp,EXPV_LEFT(left),
			   expv_int_term(INT_CONSTANT,EXPV_TYPE(right),
				     EXPV_INT_VALUE(right)+
				     EXPV_INT_VALUE((expv)EXPV_RIGHT(left)))));
	if (IS_NUMERIC_CONST_V(left) &&
	    IS_NUMERIC_CONST_V(right)) {
	    return expv_numeric_const_reduce(left, right, code, v);
	}
	break;

    case MINUS_EXPR:
	if(EXPV_IS_ZERO(left) || EXPV_IS_ZERO_FLOAT(left))
	  {
	      /* 0 - x -> unary minus */
	      code = UNARY_MINUS_EXPR;
	      left = right;
	      right = NULL;
	      break;
	  }
	if(EXPV_IS_ZERO(right) || EXPV_IS_ZERO_FLOAT(right))
	    return(left);		/* x - 0 = x */

	if(EXPV_CODE(left) == INT_CONSTANT && 
	   EXPV_CODE(right) == INT_CONSTANT)
	  return(expv_int_term(INT_CONSTANT,EXPV_TYPE(v),
			       EXPV_INT_VALUE(left)-EXPV_INT_VALUE(right)));
#ifdef REDUCE_FLOAT
	if(EXPV_CODE(left) == FLOAT_CONSTANT && 
	   EXPV_CODE(right) == FLOAT_CONSTANT)
	  return(expv_float_term(FLOAT_CONSTANT,EXPV_TYPE(v),
				 EXPV_FLOAT_VALUE(left)-EXPV_FLOAT_VALUE(right)));

#endif
	if (IS_NUMERIC_CONST_V(left) &&
	    IS_NUMERIC_CONST_V(right)) {
	    return expv_numeric_const_reduce(left, right, code, v);
	}
	break;

    case UNARY_MINUS_EXPR:
	if(EXPV_CODE(left) == INT_CONSTANT)
	  return(expv_int_term(INT_CONSTANT,EXPV_TYPE(v),
			       -EXPV_INT_VALUE(left)));
#ifdef REDUCE_FLOAT
	if(EXPV_CODE(left) == FLOAT_CONSTANT)
	  return(expv_float_term(FLOAT_CONSTANT,EXPV_TYPE(v),
				 -EXPV_FLOAT_VALUE(left)));
#endif
	if (IS_NUMERIC_CONST_V(left)) {
	    return expv_numeric_const_reduce(left, (expv)NULL, code, v);
	}
	break;

    case MOD_EXPR:
	if(EXPV_IS_ONE(right)) return(expv_constant_0); /* x%1 = 0 */
	if(EXPV_IS_MINUS_ONE(right)) return(expv_constant_0); /* x%-1 = 0 */
	if(EXPV_IS_ZERO(right)) /* x%0 = error */
	  {
	      error_at_node((expr)v, "mod by zero");
	      return(v);
	  }
	if(EXPV_CODE(left) == INT_CONSTANT && 
	   EXPV_CODE(right) == INT_CONSTANT)
	  {
	      return(expv_int_term(INT_CONSTANT,tp,
				   EXPV_INT_VALUE(left)%EXPV_INT_VALUE(right)));
	  }
	if (IS_NUMERIC_CONST_V(left) &&
	    IS_NUMERIC_CONST_V(right)) {
	    return expv_numeric_const_reduce(left, right, code, v);
	}
	break;

    case LSHIFT_EXPR:
	if(EXPV_IS_ZERO(right)) /* x<<0 = x */
	  return(left);
	if(EXPV_CODE(left) == INT_CONSTANT && 
	   EXPV_CODE(right) == INT_CONSTANT)
	  return(expv_int_term(INT_CONSTANT,tp,
			       EXPV_INT_VALUE(left)<<EXPV_INT_VALUE(right)));
	if (IS_NUMERIC_CONST_V(left) &&
	    IS_NUMERIC_CONST_V(right)) {
	    return expv_numeric_const_reduce(left, right, code, v);
	}
	break;

    case RSHIFT_EXPR:
	if(EXPV_IS_ZERO(right)) /* x>>0 = x */
	  return(left);

	if(EXPV_CODE(left) == INT_CONSTANT && 
	   EXPV_CODE(right) == INT_CONSTANT)
	  {
	      if(code == RSHIFT_EXPR)
		return(expv_int_term(INT_CONSTANT,tp,
				     EXPV_UINT_VALUE(left) >> 
				     EXPV_INT_VALUE(right)));
	      return(expv_int_term(INT_CONSTANT,tp,
				   EXPV_INT_VALUE(left)>>EXPV_INT_VALUE(right)));
	  }
	if (IS_NUMERIC_CONST_V(left) &&
	    IS_NUMERIC_CONST_V(right)) {
	    return expv_numeric_const_reduce(left, right, code, v);
	}
	break;

    case BIT_AND_EXPR:
	if(EXPV_IS_ZERO(right)||EXPV_IS_ZERO(left)) /* x & 0 = 0 */
	  return(expv_constant_0);
	if(EXPV_CODE(left) == INT_CONSTANT && 
	   EXPV_CODE(right) == INT_CONSTANT)
	  return(expv_int_term(INT_CONSTANT,tp,
			       EXPV_INT_VALUE(left)&EXPV_INT_VALUE(right)));
	if (IS_NUMERIC_CONST_V(left) &&
	    IS_NUMERIC_CONST_V(right)) {
	    return expv_numeric_const_reduce(left, right, code, v);
	}
	break;

    case BIT_OR_EXPR:
	if(EXPV_IS_MINUS_ONE(right)||EXPV_IS_MINUS_ONE(left)) 
	  return(expv_constant_m1); /* x * -1 = -1 */
	if(EXPV_CODE(left) == INT_CONSTANT && 
	   EXPV_CODE(right) == INT_CONSTANT)
	  return(expv_int_term(INT_CONSTANT,tp,
			       EXPV_INT_VALUE(left)|EXPV_INT_VALUE(right)));
	if (IS_NUMERIC_CONST_V(left) &&
	    IS_NUMERIC_CONST_V(right)) {
	    return expv_numeric_const_reduce(left, right, code, v);
	}
	break;

    case BIT_XOR_EXPR:
	if(EXPV_CODE(left) == INT_CONSTANT && 
	   EXPV_CODE(right) == INT_CONSTANT)
	  return(expv_int_term(INT_CONSTANT,tp,
			       EXPV_INT_VALUE(left)^EXPV_INT_VALUE(right)));
	if (IS_NUMERIC_CONST_V(left) &&
	    IS_NUMERIC_CONST_V(right)) {
	    return expv_numeric_const_reduce(left, right, code, v);
	}
	break;
	  
    case BIT_NOT_EXPR:
	if(EXPV_CODE(left) == INT_CONSTANT)
	  return(expv_int_term(INT_CONSTANT,tp,~EXPV_INT_VALUE(left)));
	if (IS_NUMERIC_CONST_V(left)) {
	    return expv_numeric_const_reduce(left, (expv)NULL, code, v);
	}
	break;
	  
    case LOG_EQ_EXPR:
	if(EXPV_CODE(left) == INT_CONSTANT && 
	   EXPV_CODE(right) == INT_CONSTANT)
	  return(expv_int_term(INT_CONSTANT,tp,
			       EXPV_INT_VALUE(left) == EXPV_INT_VALUE(right)));
	if(EXPV_CODE(left) == FLOAT_CONSTANT && 
	   EXPV_CODE(right) == FLOAT_CONSTANT)
	  return(expv_int_term(INT_CONSTANT,tp,
			       EXPV_FLOAT_VALUE(left)==EXPV_FLOAT_VALUE(right)));
	if (IS_NUMERIC_CONST_V(left) &&
	    IS_NUMERIC_CONST_V(right)) {
	    return expv_numeric_const_reduce(left, right, code, v);
	}
	break;

    case LOG_NEQ_EXPR:
	if(EXPV_CODE(left) == INT_CONSTANT && 
	   EXPV_CODE(right) == INT_CONSTANT)
	  return(expv_int_term(INT_CONSTANT,tp,
			       EXPV_INT_VALUE(left) != EXPV_INT_VALUE(right)));
	if(EXPV_CODE(left) == FLOAT_CONSTANT && 
	   EXPV_CODE(right) == FLOAT_CONSTANT)
	  return(expv_int_term(INT_CONSTANT,tp,
			       EXPV_FLOAT_VALUE(left)!=EXPV_FLOAT_VALUE(right)));
	if (IS_NUMERIC_CONST_V(left) &&
	    IS_NUMERIC_CONST_V(right)) {
	    return expv_numeric_const_reduce(left, right, code, v);
	}
	break;

    case LOG_GE_EXPR:
	if(EXPV_CODE(left) == INT_CONSTANT && 
	   EXPV_CODE(right) == INT_CONSTANT){
	    return(expv_int_term(INT_CONSTANT,tp,
				 EXPV_INT_VALUE(left)>=EXPV_INT_VALUE(right)));
	}
	if(EXPV_CODE(left) == FLOAT_CONSTANT && 
	   EXPV_CODE(right) == FLOAT_CONSTANT)
	  return(expv_int_term(INT_CONSTANT,tp,
			       EXPV_FLOAT_VALUE(left)>=EXPV_FLOAT_VALUE(right)));
	if (IS_NUMERIC_CONST_V(left) &&
	    IS_NUMERIC_CONST_V(right)) {
	    return expv_numeric_const_reduce(left, right, code, v);
	}
	break;

    case LOG_GT_EXPR:
	if(EXPV_CODE(left) == INT_CONSTANT && 
	   EXPV_CODE(right) == INT_CONSTANT) {
	    return(expv_int_term(INT_CONSTANT,tp,
				 EXPV_INT_VALUE(left)>EXPV_INT_VALUE(right)));
	}
	if(EXPV_CODE(left) == FLOAT_CONSTANT && 
	   EXPV_CODE(right) == FLOAT_CONSTANT)
	  return(expv_int_term(INT_CONSTANT,tp,
			       EXPV_FLOAT_VALUE(left)>EXPV_FLOAT_VALUE(right)));
	if (IS_NUMERIC_CONST_V(left) &&
	    IS_NUMERIC_CONST_V(right)) {
	    return expv_numeric_const_reduce(left, right, code, v);
	}
	break;

    case LOG_LE_EXPR:
	if(EXPV_CODE(left) == INT_CONSTANT && 
	   EXPV_CODE(right) == INT_CONSTANT){
	    return(expv_int_term(INT_CONSTANT,tp,
				 EXPV_INT_VALUE(left)<=EXPV_INT_VALUE(right)));
	}
	if(EXPV_CODE(left) == FLOAT_CONSTANT && 
	   EXPV_CODE(right) == FLOAT_CONSTANT)
	  return(expv_int_term(INT_CONSTANT,tp,
			       EXPV_FLOAT_VALUE(left)<=EXPV_FLOAT_VALUE(right)));
	if (IS_NUMERIC_CONST_V(left) &&
	    IS_NUMERIC_CONST_V(right)) {
	    return expv_numeric_const_reduce(left, right, code, v);
	}
	break;

    case LOG_LT_EXPR:
	if(EXPV_CODE(left) == INT_CONSTANT && 
	   EXPV_CODE(right) == INT_CONSTANT)
	  {
#ifdef xxx			/* problem??? */
	      return(expv_term(INT_CONSTANT,tp,NULL,
			       EXPV_UINT_VALUE(left)<
			       EXPV_UINT_VALUE(right)));
#endif
	      return(expv_int_term(INT_CONSTANT,tp,
				   EXPV_INT_VALUE(left)<EXPV_INT_VALUE(right)));
	  }
	if(EXPV_CODE(left) == FLOAT_CONSTANT &&
	   EXPV_CODE(right) == FLOAT_CONSTANT)
	  return(expv_int_term(INT_CONSTANT,tp,
			       EXPV_FLOAT_VALUE(left)<EXPV_FLOAT_VALUE(right)));
	if (IS_NUMERIC_CONST_V(left) &&
	    IS_NUMERIC_CONST_V(right)) {
	    return expv_numeric_const_reduce(left, right, code, v);
	}
	break;

    case LOG_AND_EXPR:
	if(EXPV_CODE(left) == INT_CONSTANT &&
	   EXPV_CODE(right) == INT_CONSTANT)
	  return(expv_int_term(INT_CONSTANT,tp,
			       EXPV_INT_VALUE(left) && EXPV_INT_VALUE(right)));
	if (IS_NUMERIC_CONST_V(left) &&
	    IS_NUMERIC_CONST_V(right)) {
	    return expv_numeric_const_reduce(left, right, code, v);
	}
	break;

    case LOG_OR_EXPR:
	if(EXPV_CODE(left) == INT_CONSTANT &&
	   EXPV_CODE(right) == INT_CONSTANT)
	  return(expv_int_term(INT_CONSTANT,tp,
			       EXPV_INT_VALUE(left) || EXPV_INT_VALUE(right)));
	if (IS_NUMERIC_CONST_V(left) &&
	    IS_NUMERIC_CONST_V(right)) {
	    return expv_numeric_const_reduce(left, right, code, v);
	}
	break;

    case LOG_NOT_EXPR:
	if(EXPV_CODE(left) == INT_CONSTANT)
	  return(expv_int_term(INT_CONSTANT,tp,!EXPV_INT_VALUE(left)));
	if (IS_NUMERIC_CONST_V(left)) {
	    return expv_numeric_const_reduce(left, (expv)NULL, code, v);
	}
	break;

    case FUNCTION_CALL:		/* call node, check in-line/built-in */
	if((rv = expv_inline_function(left,right)) != NULL)
	  return(rv);
    default: {}
    }
    if(code == EXPV_CODE(v) && tp == EXPV_TYPE(v) &&
       left == EXPV_LEFT(v) && right == EXPV_RIGHT(v))
      return(v);		/* no change */
    else			/* re-construct */
      return(expv_cons(code,tp,left,right));
}

#ifdef REDUCE_FLOAT
expv expv_float_reduce(expv v){ return v; }
#else
expv
expv_float_reduce(v)
     expv v;
{
    enum expr_code code;
    TYPE_DESC tp;
    expv left,right,rv;
    
    if(v == NULL) return(v);	/* error recovery */
    code = EXPV_CODE(v);
    tp = EXPV_TYPE(v);

    /* check for terminal */
    if(EXPR_CODE_IS_TERMINAL(code)) return v;
    if(EXPR_CODE_IS_LIST(code)) return v;

    /* internal node */
    left = expv_reduce(EXPV_LEFT(v));
    right = EXPV_RIGHT(v);
    if(right != NULL) right = expv_reduce(right);

    /* constant folding */
    switch(code){
    case MUL_EXPR:
	if(EXPV_CODE(left) == FLOAT_CONSTANT &&
	   EXPV_CODE(right) == FLOAT_CONSTANT)
	  return(expv_float_term(FLOAT_CONSTANT,EXPV_TYPE(v),
				 EXPV_FLOAT_VALUE(left)*EXPV_FLOAT_VALUE(right)));
	break;

    case DIV_EXPR:
	if(EXPV_CODE(left) == FLOAT_CONSTANT &&
	   EXPV_CODE(right) == FLOAT_CONSTANT)
	  return(expv_float_term(FLOAT_CONSTANT,EXPV_TYPE(v),
				 EXPV_FLOAT_VALUE(left)/EXPV_FLOAT_VALUE(right)));
	break;

    case PLUS_EXPR:
	if(EXPV_CODE(left) == FLOAT_CONSTANT && 
	   EXPV_CODE(right) == FLOAT_CONSTANT)
	  return(expv_float_term(FLOAT_CONSTANT,EXPV_TYPE(v),
				 EXPV_FLOAT_VALUE(left)+EXPV_FLOAT_VALUE(right)));
	
	break;
    case MINUS_EXPR:
	if(EXPV_CODE(left) == FLOAT_CONSTANT && 
	   EXPV_CODE(right) == FLOAT_CONSTANT)
	  return(expv_float_term(FLOAT_CONSTANT,EXPV_TYPE(v),
				 EXPV_FLOAT_VALUE(left)-EXPV_FLOAT_VALUE(right)));
	break;

    case UNARY_MINUS_EXPR:
	if(EXPV_CODE(left) == FLOAT_CONSTANT)
	  return(expv_float_term(FLOAT_CONSTANT,EXPV_TYPE(v),
				 -EXPV_FLOAT_VALUE(left)));
	break;
    }

    if(code == EXPV_CODE(v) && tp == EXPV_TYPE(v) &&
       left == EXPV_LEFT(v) && right == EXPV_RIGHT(v))
      return(v);		/* no change */
    else			/* re-construct */
      return(expv_cons(code,tp,left,right));
}
#endif


#ifdef HAS_INT64
_omInt64_t
expr_int64_value(x)
     expr x;
{
    _omInt64_t h = (_omInt64_t)(EXPR_LLINT_HIGH(x));
    _omInt64_t l = (_omInt64_t)(EXPR_LLINT_LOW(x) & INT_MASK);
    _omInt64_t i64 = (_omInt64_t)((h << INT_BITS) | (_omInt64_t)(l & INT_MASK));
#if 0
    fprintf(stderr, "debug64: val = %llx, h = %llx, l = %llx\n",
	    i64, h, l);
#endif
    return i64;
}

void
expr_int64_save(x, i64)
     expr x;
     _omInt64_t i64;
{
    int h = (int)(((int)(i64 >> INT_BITS)) & INT_MASK);
    int l = (int)(i64 & INT_MASK);

    EXPR_LLINT_HIGH(x) = h;
    EXPR_LLINT_LOW(x) = l;
#if 0
    fprintf(stderr, "debug64(sav): val = %llx, h = %x, l = %x\n",
	    i64, h, l);
#endif
}
#endif /* HAS_INT64 */    


/* 
 * convert numeric value v to type 'tp' 
 */
expv expv_reduce_conv_const(TYPE_DESC tp, expv v)
{
    if (!IS_ELEMENT_TYPE(tp)) {
	fatal("exvp_reduce_conv_const: bad type");
	return NULL;
    }

    if (tp == NULL) {
	return NULL;
    }

    if (EXPV_TYPE(v) == NULL) {
	return NULL;
    }

    if (TYPE_BASIC_TYPE(tp) == TYPE_BASIC_TYPE(EXPV_TYPE(v))) {
	return v;
    }

#if 0
    fprintf(stderr, "debug: tp='%s' v='%s'\n",
	    basic_type_name(TYPE_BASIC_TYPE(tp)),
	    basic_type_name(TYPE_BASIC_TYPE(EXPV_TYPE(v))));
    expr_print(v, stderr);
#endif

    switch (TYPE_BASIC_TYPE(tp)) {

#ifdef ENABLE_QREAL
	case TYPE_QREAL: {
	    return expv_reduce_convTo_qreal_const(tp, v);
	}
#endif /* ENABLE_QREAL */

	case TYPE_REAL: {
	    float f;
	    switch (TYPE_BASIC_TYPE(EXPV_TYPE(v))) {
#ifdef HAS_INT64
		case TYPE_LONGLONG: {
		    _omInt64_t i64 = EXPV_INT64_VALUE(v);
		    f = (float)i64;
		    break;
		}
#endif /* HAS_INT64 */
		case TYPE_INT: {
		    f = (float)EXPV_INT_VALUE(v);
		    break;
		}
		case TYPE_DREAL:
		case TYPE_REAL: {
		    f = (float)EXPV_FLOAT_VALUE(v);
		    break;
		}
		case TYPE_DCOMPLEX: {
		    expv re = expv_complex_real(type_DREAL, v);
		    f = (float)EXPV_FLOAT_VALUE(re);
		    break;
		}
		case TYPE_COMPLEX: {
		    expv re = expv_complex_real(type_REAL, v);
		    f = (float)EXPV_FLOAT_VALUE(re);
		    break;
		}
#ifdef ENABLE_QREAL
		case TYPE_QREAL: {
		    return expv_reduce_convFrom_qreal_const(tp, v);
		}
#endif /* ENABLE_QREAL */

		default: {
		    fatal("expv_reduce_conv_const: not a numeric constant.");
		    return NULL;
		}
	    }
	    return expv_float_term(FLOAT_CONSTANT, tp, (double)f);
	}

	case TYPE_DREAL: {
	    double d;
	    switch (TYPE_BASIC_TYPE(EXPV_TYPE(v))) {
#ifdef HAS_INT64
		case TYPE_LONGLONG: {
		    _omInt64_t i64 = EXPV_INT_VALUE(v);
		    d = (double)i64;
		    break;
		}
#endif /* HAS_INT64 */
		case TYPE_INT: {
		    d = (double)EXPV_INT_VALUE(v);
		    break;
		}
		case TYPE_DREAL:
		case TYPE_REAL: {
		    d = (double)EXPV_FLOAT_VALUE(v);
		    break;
		}
		case TYPE_DCOMPLEX: {
		    expv re = expv_complex_real(type_DREAL, v);
		    d = (double)EXPV_FLOAT_VALUE(re);
		    break;
		}
		case TYPE_COMPLEX: {
		    expv re = expv_complex_real(type_REAL, v);
		    d = (double)EXPV_FLOAT_VALUE(re);
		    break;
		}
#ifdef ENABLE_QREAL
		case TYPE_QREAL: {
		    return expv_reduce_convFrom_qreal_const(tp, v);
		}
#endif /* ENABLE_QREAL */
		default: {
		    fatal("expv_reduce_conv_const: not a numeric constant.");
		    return NULL;
		}
	    }
	    return expv_float_term(FLOAT_CONSTANT, tp, (double)d);
	}

	case TYPE_CHAR:
	case TYPE_SHORT:
	case TYPE_INT: {
	    int i;
	    switch (TYPE_BASIC_TYPE(EXPV_TYPE(v))) {
#ifdef HAS_INT64
		case TYPE_LONGLONG: {
		    _omInt64_t i64 = EXPV_INT64_VALUE(v);
		    i = (int)i64;
		    break;
		}
#endif /* HAS_INT64 */
		case TYPE_INT: {
		    i = (int)EXPV_INT_VALUE(v);
		    break;
		}
		case TYPE_DREAL:
		case TYPE_REAL: {
		    i = (int)EXPV_FLOAT_VALUE(v);
		    break;
		}
		case TYPE_DCOMPLEX: {
		    expv re = expv_complex_real(type_DREAL, v);
		    i = (int)EXPV_FLOAT_VALUE(re);
		    break;
		}
		case TYPE_COMPLEX: {
		    expv re = expv_complex_real(type_REAL, v);
		    i = (int)EXPV_FLOAT_VALUE(re);
		    break;
		}
#ifdef ENABLE_QREAL
		case TYPE_QREAL: {
		    return expv_reduce_convFrom_qreal_const(tp, v);
		}
#endif /* ENABLE_QREAL */
		default: {
		    fatal("expv_reduce_conv_const: not a numeric constant.");
		    return NULL;
		}
	    }
	    if (TYPE_BASIC_TYPE(tp) == TYPE_CHAR) {
		unsigned char c = (unsigned char)i;
		return expv_int_term(INT_CONSTANT, tp, (int)c);
	    } else if (TYPE_BASIC_TYPE(tp) == TYPE_SHORT) {
		unsigned short s = (unsigned short)i;
		return expv_int_term(INT_CONSTANT, tp, (int)s);
	    } else {
		return expv_int_term(INT_CONSTANT, tp, i);
	    }
	}

#ifdef HAS_INT64
	case TYPE_LONGLONG: {
	    _omInt64_t i64 = (_omInt64_t)0;
	    expv ret = expv_longlong_term(LONGLONG_CONSTANT, tp, 0, 0);
	    switch (TYPE_BASIC_TYPE(EXPV_TYPE(v))) {
		case TYPE_LONGLONG: {
		    i64 = (_omInt64_t)EXPV_INT64_VALUE(v);
		    break;
		}
		case TYPE_INT: {
		    i64 = (_omInt64_t)EXPV_INT_VALUE(v);
		    break;
		}
		case TYPE_DREAL:
		case TYPE_REAL: {
		    i64 = (_omInt64_t)EXPV_FLOAT_VALUE(v);
		    break;
		}
		case TYPE_DCOMPLEX: {
		    expv re = expv_complex_real(type_DREAL, v);
		    i64 = (_omInt64_t)EXPV_FLOAT_VALUE(re);
		    break;
		}
		case TYPE_COMPLEX: {
		    expv re = expv_complex_real(type_REAL, v);
		    i64 = (_omInt64_t)EXPV_FLOAT_VALUE(re);
		    break;
		}
#ifdef ENABLE_QREAL
		case TYPE_QREAL: {
		    return expv_reduce_convFrom_qreal_const(tp, v);
		}
#endif /* ENABLE_QREAL */
		default: {
		    fatal("expv_reduce_conv_const: not a numeric constant.");
		    return NULL;
		}
	    }
	    EXPV_INT64_SAVE_VALUE(ret, i64);
	    return ret;
	}
#endif /* HAS_INT64 */

	case TYPE_COMPLEX: {
	    float f;
	    switch (TYPE_BASIC_TYPE(EXPV_TYPE(v))) {
#ifdef HAS_INT64
		case TYPE_LONGLONG: {
		    _omInt64_t i64 = EXPV_INT64_VALUE(v);
		    f = (float)i64;
		    break;
		}
#endif /* HAS_INT64 */
		case TYPE_INT: {
		    f = (float)EXPV_INT_VALUE(v);
		    break;
		}
		case TYPE_DREAL:
		case TYPE_REAL: {
		    f = (float)EXPV_FLOAT_VALUE(v);
		    break;
		}
		case TYPE_COMPLEX:
		case TYPE_DCOMPLEX: {
		    expv re = expv_complex_real(type_DREAL, v);
		    expv im = expv_complex_img(type_DREAL, v);
		    return expv_z_cons(re, im, TRUE);
		}
#ifdef ENABLE_QREAL
		case TYPE_QREAL: {
		    return expv_reduce_convFrom_qreal_const(tp, v);
		}
#endif /* ENABLE_QREAL */
		default: {
		    fatal("expv_reduce_conv_const: not a numeric constant.");
		    return NULL;
		}
	    }
	    return expv_c_cons(expv_float_term(FLOAT_CONSTANT, type_REAL, (double)f),
			       expv_float_0,
			       TRUE);
	}
    
	case TYPE_DCOMPLEX: {
	    double d;
	    switch (TYPE_BASIC_TYPE(EXPV_TYPE(v))) {
#ifdef HAS_INT64
		case TYPE_LONGLONG: {
		    _omInt64_t i64 = EXPV_INT64_VALUE(v);
		    d = (double)i64;
		    break;
		}
#endif /* HAS_INT64 */
		case TYPE_INT: {
		    d = (double)EXPV_INT_VALUE(v);
		    break;
		}
		case TYPE_DREAL:
		case TYPE_REAL: {
		    d = (double)EXPV_FLOAT_VALUE(v);
		    break;
		}
		case TYPE_COMPLEX:
		case TYPE_DCOMPLEX: {
		    expv re = expv_complex_real(type_REAL, v);
		    expv im = expv_complex_img(type_REAL, v);
		    return expv_c_cons(re, im, TRUE);
		}
#ifdef ENABLE_QREAL
		case TYPE_QREAL: {
		    return expv_reduce_convFrom_qreal_const(tp, v);
		}
#endif /* ENABLE_QREAL */
		default: {
		    fatal("expv_reduce_conv_const: not a numeric constant.");
		    return NULL;
		}
	    }
	    return expv_z_cons(expv_float_term(FLOAT_CONSTANT, type_DREAL, d),
			       expv_double_0,
			       TRUE);
	}

	default: {
	    fatal("expv_reduce_conv_const: bad arithmetic type");
	}
    }
    return NULL;
}


#ifdef ENABLE_QREAL
expv
expv_reduce_convTo_qreal_const(tp, v)
     TYPE_DESC tp;
     expv v;
{
    mpf_t val;
    expv ret = NULL;
    char *str = NULL;
    
    if (!(IS_NUMERIC_CONST_V(v))) {
	error("about to convert non numeric to quad real.");
	return NULL;
    }

    if (EXPV_CODE(v) == QREAL_CONSTANT) {
	return v;
    }

    mpf_init(val);

    if (expv_const_convert_to_mpf(v, val) == FALSE) {
	fatal("expv_reduce_convTo_qreal_const: failed to convert to mpf_t??");
	mpf_clear(val);
	return NULL;
    }

    str = mpf2str(val);
    ret = expv_qreal_term(QREAL_CONSTANT, type_QREAL,
			  expv_str_term(STRING_CONSTANT,
					type_char(strlen(str)),
					str));
    mpf_clear(val);
    return ret;
}


expv
expv_reduce_convFrom_qreal_const(tp, v)
     TYPE_DESC tp;
     expv v;
{
    mpf_t val;
    double dV = 0.0;
    char *str = NULL;
    expv ret = NULL;

    if (TYPE_BASIC_TYPE(EXPV_TYPE(v)) != TYPE_QREAL &&
	EXPV_CODE(v) != QREAL_CONSTANT &&
	EXPV_CODE(v) != F_QREAL_CONSTANT) {
	fatal("expv_reduce_convFrom_qreal_const: not a quad-real constant.");
	return NULL;
    }
    
    if (TYPE_BASIC_TYPE(tp) == TYPE_QREAL &&
	EXPV_CODE(v) == QREAL_CONSTANT) {
	return v;
    }

    if (EXPV_CODE(v) == QREAL_CONSTANT) {
	str = EXPV_STR(((expv)EXPV_QREAL_CONST(v)));
    } else {
	str = (char *)EXPR_QREAL_CONST(v);
    }

    if (mpf_init_set_str(val, str, 10) < 0) {
	fatal("expv_reduce_convFrom_qreal_const: failed to assign a val to mpf_t??");
	mpf_clear(val);
	return NULL;
    }

    dV = mpf_get_d(val);
    mpf_clear(val);
    
    switch (TYPE_BASIC_TYPE(tp)) {

#ifdef HAS_INT64
	case TYPE_LONGLONG: {
	    _omInt64_t i64 = (_omInt64_t)dV;
	    ret = expv_longlong_term(LONGLONG_CONSTANT, tp, 0, 0);
	    EXPV_INT64_SAVE_VALUE(ret, i64);
	    break;
	}
#endif /* HAS_INT64 */

	case TYPE_SHORT:
	case TYPE_INT:
	case TYPE_LOGICAL: {
	    int i = (int)dV;
	    ret = expv_int_term(INT_CONSTANT, tp, i);
	    break;
	}

	case TYPE_REAL:
	case TYPE_DREAL: {
	    ret = expv_float_term(FLOAT_CONSTANT, tp, dV);
	    break;
	}

	case TYPE_COMPLEX: {
	    ret = expv_cons(COMPLEX_CONSTANT, tp,
			    expv_float_term(FLOAT_CONSTANT, type_REAL, dV),
			    expv_float_0);
	    break;
	}
	case TYPE_DCOMPLEX: {
	    ret = expv_cons(COMPLEX_CONSTANT, tp,
			    expv_float_term(FLOAT_CONSTANT, type_DREAL, dV),
			    expv_double_0);
	    break;
	}

	default: {
	    fatal("expv_reduce_convFrom_qreal_const: unknown type.");
	    break;
	}
    }

    return ret;
}
#endif /* ENABLE_QREAL */


expv expv_inline_function(left,right)
     expv left,right;
{
#ifdef not
    TYPE_DESC tp;
    char *name;
    int d,incr;
    
    if(EXPV_CODE(left) != LABEL_CONSTANT) return(NULL);
    name = EXPV_NAME(left);
    if(strcmp(name,"__builtin_va_arg_incr") == 0)
      {
	  if(EXPV_CODE(right) == ARGLIST ||
	     !IS_POINTER(EXPV_TYPE(right)))
	    {
		error("bad usage, __builtin_va_arg_incr");
		return(NULL);
	    }
	  tp = EXPV_TYPE(right);
	  d = TYPE_SIZE(TYPE_REF(tp));
	  incr = ROUND(d,PARAM_ALIGN);
	  /* __builtin_va_arg_incr(ap) => 
	     (MINUS_EXPR (INCR ap ARG_SIZE) TYPE_SIZE) */
	  if(d == incr)
	    return(expv_cons(POST_INCR_EXPR,tp,right,
			     expv_term(INT_CONSTANT,int_type,NULL,d)));
	  else
	    return(expv_cons(MINUS_EXPR,tp,
			     expv_cons(ASG_PLUS_EXPR,tp,right,
				       expv_term(INT_CONSTANT,int_type,NULL,
						 incr)),
			     expv_term(INT_CONSTANT,int_type,NULL,d)));
      }
    else 
#endif
      return(NULL);
}

/* borrowed from libf2c */
double power_di(double x,int n)
{
    double pow;
    unsigned long u;

    pow = 1;
    if(n != 0){
        if(n < 0){
	    n = -n;
	    x = 1/x;
	}
        for(u = n; ; ){
	    if(u & 01)
		pow *= x;
	    if(u >>= 1)
		x *= x;
	    else
		break;
	}
    }
    return(pow);
}

int power_ii(int x, int n)
{
    int pow;
    unsigned long u;

    if (n <= 0) {
	if (n == 0 || x == 1)
	    return 1;
	if (x != -1)
	    return x == 0 ? 1/x : 0;
	n = -n;
    }
    u = n;
    for(pow = 1; ; )
    {
	if(u & 01)
	    pow *= x;
	if(u >>= 1)
	    x *= x;
	else
	    break;
    }
    return(pow);
}

static void z_div(doublecomplex *c, doublecomplex *a, doublecomplex *b)
{
    double ratio, den;
    double abr, abi, cr;

    if ((abr = b->r) < 0.) {
	abr = - abr;
    }
    if ((abi = b->i) < 0.) {
	abi = - abi;
    }
    if ( abr <= abi ) {
	if (abi == 0) {
	    error("complex division by zero");
	}
	ratio = b->r / b->i;
	den = b->i * (1 + ratio*ratio);
	cr = (a->r*ratio + a->i) / den;
	c->i = (a->i*ratio - a->r) / den;
    } else {
	ratio = b->i / b->r ;
	den = b->r * (1 + ratio*ratio);
	cr = (a->r + a->i*ratio) / den;
	c->i = (a->i - a->r*ratio) / den;
    }
    c->r = cr;
}

void power_zi(doublecomplex *p, doublecomplex *a, int b)
{
    int n;
    unsigned long u;
    double t;
    doublecomplex q, x;
    static doublecomplex one = {1.0, 0.0};

    n = b;
    q.r = 1;
    q.i = 0;

    if (n == 0) {
	goto done;
    }
    if (n < 0) {
	n = -n;
	z_div(&x, &one, a);
    } else {
	x.r = a->r;
	x.i = a->i;
    }

    for (u = n; ; ) {
	if (u & 01) {
	    t = q.r * x.r - q.i * x.i;
	    q.i = q.r * x.i + q.i * x.r;
	    q.r = t;
	}
	if (u >>= 1) {
	    t = x.r * x.r - x.i * x.i;
	    x.i = 2 * x.r * x.i;
	    x.r = t;
	} else {
	    break;
	}
    }
    done:
    p->i = q.i;
    p->r = q.r;
}

static double f__cabs(double real, double imag)
{
    double temp;

    if (real < 0) {
	real = -real;
    }
    if (imag < 0) {
	imag = -imag;
    }
    if (imag > real) {
	temp = real;
	real = imag;
	imag = temp;
    }
    if ((real + imag) == real) {
	return real;
    }

    temp = imag/real;
    temp = real*sqrt(1.0 + temp*temp);  /*overflow!!*/
    return temp;
}

void power_zz(doublecomplex *r, doublecomplex *a, doublecomplex *b)
{
    double logr, logi, x, y;

    logr = log( f__cabs(a->r, a->i) );
    logi = atan2(a->i, a->r);

    x = exp( logr * b->r - logi * b->i );
    y = logr * b->i + logi * b->r;
    
    r->r = x * cos(y);
    r->i = x * sin(y);
}


static int complexConstFolding _ANSI_ARGS_((expv v, TYPE_DESC tp, doublecomplex *zPtr));

static int
complexConstFolding(v, tp, zPtr)
     expv v;
     TYPE_DESC tp;
     doublecomplex *zPtr;
{
    TYPE_DESC fTp = (TYPE_BASIC_TYPE(tp) == TYPE_DCOMPLEX) ? type_DREAL : type_REAL;

    /* terminal */
    if (EXPR_CODE_IS_TERMINAL_OR_CONST(EXPV_CODE(v))) {
	expv rV, iV;

	v = compile_terminal_node(v);
	if (TYPE_BASIC_TYPE(EXPV_TYPE(v)) != TYPE_BASIC_TYPE(tp)) {
	    v = expv_reduce_conv_const(tp, v);
	}
	rV = expv_complex_real(fTp, v);
	iV = expv_complex_img(fTp, v);

	if (expr_is_constant(rV) == FALSE) {
	    error("real part is not a constant");
	    return FALSE;
	} else {
	    rV = expr_constant_value(rV, FALSE);
	    if (rV == NULL) {
		error("can't get constant value of real part");
		return FALSE;
	    }
	}
	if (expr_is_constant(iV) == FALSE) {
	    error("imaginary part is not a constant");
	    return FALSE;
	} else {
	    iV = expr_constant_value(iV, FALSE);
	    if (iV == NULL) {
		error("can't get constant value of imaginary part");
		return FALSE;
	    }
	}

	if (TYPE_BASIC_TYPE(EXPV_TYPE(rV)) != TYPE_BASIC_TYPE(fTp)) {
	    rV = expv_reduce_conv_const(fTp, rV);
	}
	if (TYPE_BASIC_TYPE(EXPV_TYPE(iV)) != TYPE_BASIC_TYPE(fTp)) {
	    iV = expv_reduce_conv_const(fTp, iV);
	}

	zPtr->r = (double)EXPV_FLOAT_VALUE(rV);
	zPtr->i = (double)EXPV_FLOAT_VALUE(iV);

	return TRUE;
    }

    switch (EXPV_CODE(v)) {

	case F_UNARY_MINUS_EXPR:
	case UNARY_MINUS_EXPR: {
	    doublecomplex zVal;
	    if (complexConstFolding(EXPR_ARG1(v), tp, &zVal) == FALSE) {
		return FALSE;
	    }
	    zPtr->r = -(zVal.r);
	    zPtr->i = -(zVal.i);
	    return TRUE;
	}

	case F_PLUS_EXPR:
	case PLUS_EXPR: {
	    doublecomplex z0, z1;

	    if (complexConstFolding(EXPR_ARG1(v), tp, &z0) == FALSE) {
		return FALSE;
	    }
	    if (complexConstFolding(EXPR_ARG2(v), tp, &z1) == FALSE) {
		return FALSE;
	    }

	    zPtr->r = z0.r + z1.r;
	    zPtr->i = z0.i + z1.i;
	    return TRUE;
	}

	case F_MINUS_EXPR:
	case MINUS_EXPR: {
	    doublecomplex z0, z1;

	    if (complexConstFolding(EXPR_ARG1(v), tp, &z0) == FALSE) {
		return FALSE;
	    }
	    if (complexConstFolding(EXPR_ARG2(v), tp, &z1) == FALSE) {
		return FALSE;
	    }

	    zPtr->r = z0.r - z1.r;
	    zPtr->i = z0.i - z1.i;
	    return TRUE;
	}

	case F_MUL_EXPR:
	case MUL_EXPR: {
	    doublecomplex z0, z1;

	    if (complexConstFolding(EXPR_ARG1(v), tp, &z0) == FALSE) {
		return FALSE;
	    }
	    if (complexConstFolding(EXPR_ARG2(v), tp, &z1) == FALSE) {
		return FALSE;
	    }
	    /*
	     * R = left.r * right.r - left.i * right.i
	     * I = left.r * right.i + left.i * right.r
	     */
	    zPtr->r = z0.r * z1.r - z0.i * z1.i;
	    zPtr->i = z0.r * z1.i + z0.i * z1.r;
	    return TRUE;
	}

	case F_DIV_EXPR:
	case DIV_EXPR: {
	    double norm;
	    doublecomplex z0, z1;

	    if (complexConstFolding(EXPR_ARG1(v), tp, &z0) == FALSE) {
		return FALSE;
	    }
	    if (complexConstFolding(EXPR_ARG2(v), tp, &z1) == FALSE) {
		return FALSE;
	    }
	    /*
	     * norm = right.r ^ 2 + right.i ^ 2
	     * R = (left.r * right.r + left.i * right.i) / norm
	     * I = (left.i * right.r - left.r * right.i) / norm
	     */
	    norm = z1.r * z1.r + z1.i * z1.i;
	    zPtr->r = (z0.r * z1.r + z0.i * z1.i) / norm;
	    zPtr->i = (z0.i * z1.r - z0.r * z1.i) / norm;
	    return TRUE;
	}

	case F_POWER_EXPR: {
	    doublecomplex z0, z1, zz;

	    if (complexConstFolding(EXPR_ARG1(v), tp, &z0) == FALSE) {
		return FALSE;
	    }
	    if (complexConstFolding(EXPR_ARG2(v), tp, &z1) == FALSE) {
		return FALSE;
	    }

	    power_zz(&zz, &z0, &z1);
	    zPtr->r = zz.r;
	    zPtr->i = zz.i;
	    return TRUE;
	}
	    
	default: {
	    error("illegal complex operation.");
	    return FALSE;
	}

    }

    return FALSE;
}


expv
expv_complex_const_reduce(v, tp)
     expv v;
     TYPE_DESC tp;
{
    doublecomplex z;
    expv vI, vR;

    if (complexConstFolding(v, tp, &z) == FALSE) {
	return NULL;
    }

    if (TYPE_BASIC_TYPE(tp) == TYPE_COMPLEX) {
	vR = expv_float_term(FLOAT_CONSTANT, type_REAL, z.r);
	vI = expv_float_term(FLOAT_CONSTANT, type_REAL, z.i);
    } else {
	vR = expv_float_term(FLOAT_CONSTANT, type_DREAL, z.r);
	vI = expv_float_term(FLOAT_CONSTANT, type_DREAL, z.i);
    }
    return expv_cons(COMPLEX_CONSTANT, tp, vR, vI);
}
