"Fossies" - the Fresh Open Source Software Archive

Member "mathmod-branches-r508-trunk/fparser/fparser.cc" (8 Mar 2021, 150833 Bytes) of package /linux/misc/mathmod-11.0-source.zip:


As a special service "Fossies" has tried to format the requested source page into HTML format using (guessed) C and C++ source code syntax highlighting (style: standard) with prefixed line numbers and code folding option. Alternatively you can here view or download the uninterpreted source code file. For more information about "fparser.cc" see the Fossies "Dox" file reference documentation and the latest Fossies "Diffs" side-by-side code changes report: 10.1_vs_11.0.

    1 /***************************************************************************\
    2 |* Function Parser for C++ v4.5.2                                          *|
    3 |*-------------------------------------------------------------------------*|
    4 |* Copyright: Juha Nieminen, Joel Yliluoma                                 *|
    5 |*                                                                         *|
    6 |* This library is distributed under the terms of the                      *|
    7 |* GNU Lesser General Public License version 3.                            *|
    8 |* (See lgpl.txt and gpl.txt for the license text.)                        *|
    9 \***************************************************************************/
   10 
   11 #include "fpconfig.hh"
   12 #include "fparser.hh"
   13 #include <set>
   14 #include <cstdlib>
   15 #include <cstring>
   16 #include <cctype>
   17 #include <cmath>
   18 #include <cassert>
   19 #include <limits>
   20 #include "extrasrc/fptypes.hh"
   21 #include "extrasrc/fpaux.hh"
   22 using namespace FUNCTIONPARSERTYPES;
   23 
   24 #ifdef FP_USE_THREAD_SAFE_EVAL_WITH_ALLOCA
   25 #ifndef FP_USE_THREAD_SAFE_EVAL
   26 #define FP_USE_THREAD_SAFE_EVAL
   27 #endif
   28 #endif
   29 
   30 #ifdef __GNUC__
   31 # define likely(x)       __builtin_expect(!!(x), 1)
   32 # define unlikely(x)     __builtin_expect(!!(x), 0)
   33 #else
   34 # define likely(x)   (x)
   35 # define unlikely(x) (x)
   36 #endif
   37 
   38 //=========================================================================
   39 // Opcode analysis functions
   40 //=========================================================================
   41 // These functions are used by the Parse() bytecode optimizer (mostly from
   42 // code in fp_opcode_add.inc).
   43 
   44 bool FUNCTIONPARSERTYPES::IsLogicalOpcode(unsigned op)
   45 {
   46     switch(op)
   47     {
   48       case cAnd: case cAbsAnd:
   49       case cOr:  case cAbsOr:
   50       case cNot: case cAbsNot:
   51       case cNotNot: case cAbsNotNot:
   52       case cEqual: case cNEqual:
   53       case cLess: case cLessOrEq:
   54       case cGreater: case cGreaterOrEq:
   55           return true;
   56       default: break;
   57     }
   58     return false;
   59 }
   60 
   61 bool FUNCTIONPARSERTYPES::IsComparisonOpcode(unsigned op)
   62 {
   63     switch(op)
   64     {
   65       case cEqual: case cNEqual:
   66       case cLess: case cLessOrEq:
   67       case cGreater: case cGreaterOrEq:
   68           return true;
   69       default: break;
   70     }
   71     return false;
   72 }
   73 
   74 unsigned FUNCTIONPARSERTYPES::OppositeComparisonOpcode(unsigned op)
   75 {
   76     switch(op)
   77     {
   78       case cLess: return cGreater;
   79       case cGreater: return cLess;
   80       case cLessOrEq: return cGreaterOrEq;
   81       case cGreaterOrEq: return cLessOrEq;
   82     }
   83     return op;
   84 }
   85 
   86 bool FUNCTIONPARSERTYPES::IsNeverNegativeValueOpcode(unsigned op)
   87 {
   88     switch(op)
   89     {
   90       case cAnd: case cAbsAnd:
   91       case cOr:  case cAbsOr:
   92       case cNot: case cAbsNot:
   93       case cNotNot: case cAbsNotNot:
   94       case cEqual: case cNEqual:
   95       case cLess: case cLessOrEq:
   96       case cGreater: case cGreaterOrEq:
   97       case cSqrt: case cRSqrt: case cSqr:
   98       case cHypot:
   99       case cAbs:
  100       case cAcos: case cCosh:
  101           return true;
  102       default: break;
  103     }
  104     return false;
  105 }
  106 
  107 bool FUNCTIONPARSERTYPES::IsAlwaysIntegerOpcode(unsigned op)
  108 {
  109     switch(op)
  110     {
  111       case cAnd: case cAbsAnd:
  112       case cOr:  case cAbsOr:
  113       case cNot: case cAbsNot:
  114       case cNotNot: case cAbsNotNot:
  115       case cEqual: case cNEqual:
  116       case cLess: case cLessOrEq:
  117       case cGreater: case cGreaterOrEq:
  118       case cInt: case cFloor: case cCeil: case cTrunc:
  119           return true;
  120       default: break;
  121     }
  122     return false;
  123 }
  124 
  125 bool FUNCTIONPARSERTYPES::IsUnaryOpcode(unsigned op)
  126 {
  127     switch(op)
  128     {
  129       case cInv: case cNeg:
  130       case cNot: case cAbsNot:
  131       case cNotNot: case cAbsNotNot:
  132       case cSqr: case cRSqrt:
  133       case cDeg: case cRad:
  134           return true;
  135     }
  136     return (op < FUNC_AMOUNT && Functions[op].params == 1);
  137 }
  138 
  139 bool FUNCTIONPARSERTYPES::IsBinaryOpcode(unsigned op)
  140 {
  141     switch(op)
  142     {
  143       case cAdd: case cSub: case cRSub:
  144       case cMul: case cDiv: case cRDiv:
  145       case cMod:
  146       case cEqual: case cNEqual: case cLess:
  147       case cLessOrEq: case cGreater: case cGreaterOrEq:
  148       case cAnd: case cAbsAnd:
  149       case cOr: case cAbsOr:
  150           return true;
  151     }
  152     return (op < FUNC_AMOUNT && Functions[op].params == 2);
  153 }
  154 
  155 bool FUNCTIONPARSERTYPES::IsVarOpcode(unsigned op)
  156 {
  157     // See comment in declaration of FP_ParamGuardMask
  158     return int(op) >= VarBegin;
  159 }
  160 
  161 bool FUNCTIONPARSERTYPES::IsCommutativeOrParamSwappableBinaryOpcode(unsigned op)
  162 {
  163     switch(op)
  164     {
  165       case cAdd:
  166       case cMul:
  167       case cEqual: case cNEqual:
  168       case cAnd: case cAbsAnd:
  169       case cOr: case cAbsOr:
  170       case cMin: case cMax: case cHypot:
  171           return true;
  172       case cDiv: case cSub: case cRDiv: case cRSub:
  173           return true;
  174       case cLess: case cGreater:
  175       case cLessOrEq: case cGreaterOrEq: return true;
  176     }
  177     return false;
  178 }
  179 
  180 unsigned FUNCTIONPARSERTYPES::GetParamSwappedBinaryOpcode(unsigned op)
  181 {
  182     switch(op)
  183     {
  184       case cAdd:
  185       case cMul:
  186       case cEqual: case cNEqual:
  187       case cAnd: case cAbsAnd:
  188       case cOr: case cAbsOr:
  189       case cMin: case cMax: case cHypot:
  190           return op;
  191       case cDiv: return cRDiv;
  192       case cSub: return cRSub;
  193       case cRDiv: return cDiv;
  194       case cRSub: return cSub;
  195       case cLess: return cGreater;
  196       case cGreater: return cLess;
  197       case cLessOrEq: return cGreaterOrEq;
  198       case cGreaterOrEq: return cLessOrEq;
  199     }
  200     return op; // Error
  201 }
  202 
  203 template<bool ComplexType>
  204 bool FUNCTIONPARSERTYPES::HasInvalidRangesOpcode(unsigned op)
  205 {
  206     // Returns true, if the given opcode has a range of
  207     // input values that gives an error.
  208     if(ComplexType)
  209     {
  210         // COMPLEX:
  211         switch(op)
  212         {
  213           case cAtan:  // allowed range: x != +-1i
  214           case cAtanh: // allowed range: x != +-1
  215           //case cCot: // allowed range: tan(x) != 0
  216           //case cCsc: // allowed range: sin(x) != 0
  217           case cLog:   // allowed range: x != 0
  218           case cLog2:  // allowed range: x != 0
  219           case cLog10: // allowed range: x != 0
  220     #ifdef FP_SUPPORT_OPTIMIZER
  221           case cLog2by:// allowed range: x != 0
  222     #endif
  223           //case cPow: // allowed when: x != 0 or y != 0
  224           //case cSec: // allowed range: cos(x) != 0
  225           //case cTan:   // allowed range: cos(x) != 0  --> x != +-(pi/2)
  226           //case cTanh:  // allowed range: log(x) != -1 --> x != +-(pi/2)i
  227           case cRSqrt: // allowed range: x != 0
  228           //case cDiv: // allowed range: y != 0
  229           //case cRDiv: // allowed range: x != 0
  230           //case cInv: // allowed range: x != 0
  231           return true;
  232         }
  233     }
  234     else
  235     {
  236         // REAL:
  237         switch(op)
  238         {
  239           case cAcos: // allowed range: |x| <= 1
  240           case cAsin: // allowed range: |x| <= 1
  241           case cAcosh: // allowed range: x >= 1
  242           case cAtanh: // allowed range: |x| < 1
  243           //case cCot: // allowed range: tan(x) != 0
  244           //case cCsc: // allowed range: sin(x) != 0
  245           case cLog:   // allowed range: x > 0
  246           case cLog2:  // allowed range: x > 0
  247           case cLog10: // allowed range: x > 0
  248     #ifdef FP_SUPPORT_OPTIMIZER
  249           case cLog2by:// allowed range: x > 0
  250     #endif
  251           //case cPow: // allowed when: x > 0 or (x = 0 and y != 0) or (x<0)
  252                        // Technically, when (x<0 and y is not integer),
  253                        // it is not allowed, but we allow it anyway
  254                        // in order to make nontrivial roots work.
  255           //case cSec: // allowed range: cos(x) != 0
  256           case cSqrt: // allowed range: x >= 0
  257           case cRSqrt: // allowed range: x > 0
  258           //case cTan:   // allowed range: cos(x) != 0 --> x != +-(pi/2)
  259           //case cDiv: // allowed range: y != 0
  260           //case cRDiv: // allowed range: x != 0
  261           //case cInv: // allowed range: x != 0
  262           return true;
  263         }
  264     }
  265     return false;
  266 }
  267 
  268 template bool FUNCTIONPARSERTYPES::HasInvalidRangesOpcode<false>(unsigned op);
  269 template bool FUNCTIONPARSERTYPES::HasInvalidRangesOpcode<true>(unsigned op);
  270 
  271 
  272 #if(0) // Implementation moved to fpaux.hh due to linker problems
  273 //=========================================================================
  274 // Mathematical template functions
  275 //=========================================================================
  276 /* fp_pow() is a wrapper for std::pow()
  277  * that produces an identical value for
  278  * exp(1) ^ 2.0  (0x4000000000000000)
  279  * as exp(2.0)   (0x4000000000000000)
  280  * - std::pow() on x86_64
  281  * produces 2.0  (0x3FFFFFFFFFFFFFFF) instead!
  282  * See comments below for other special traits.
  283  */
  284 namespace
  285 {
  286     template<typename ValueT>
  287     inline ValueT fp_pow_with_exp_log(const ValueT& x, const ValueT& y)
  288     {
  289         // Exponentiation using exp(log(x)*y).
  290         // See http://en.wikipedia.org/wiki/Exponentiation#Real_powers
  291         // Requirements: x > 0.
  292         return fp_exp(fp_log(x) * y);
  293     }
  294 
  295     template<typename ValueT>
  296     inline ValueT fp_powi(ValueT x, unsigned long y)
  297     {
  298         // Fast binary exponentiation algorithm
  299         // See http://en.wikipedia.org/wiki/Exponentiation_by_squaring
  300         // Requirements: y is non-negative integer.
  301         ValueT result(1);
  302         while(y != 0)
  303         {
  304             if(y & 1) { result *= x; y -= 1; }
  305             else      { x *= x;      y /= 2; }
  306         }
  307         return result;
  308     }
  309 }
  310 
  311 template<typename ValueT>
  312 ValueT FUNCTIONPARSERTYPES::fp_pow(const ValueT& x, const ValueT& y)
  313 {
  314     if(x == ValueT(1)) return ValueT(1);
  315     // y is now zero or positive
  316     if(isLongInteger(y))
  317     {
  318         // Use fast binary exponentiation algorithm
  319         if(y >= ValueT(0))
  320             return fp_powi(x,              makeLongInteger(y));
  321         else
  322             return ValueT(1) / fp_powi(x, -makeLongInteger(y));
  323     }
  324     if(y >= ValueT(0))
  325     {
  326         // y is now positive. Calculate using exp(log(x)*y).
  327         if(x > ValueT(0)) return fp_pow_with_exp_log(x, y);
  328         if(x == ValueT(0)) return ValueT(0);
  329         // At this point, y > 0.0 and x is known to be < 0.0,
  330         // because positive and zero cases are already handled.
  331         if(!isInteger(y*ValueT(16)))
  332             return -fp_pow_with_exp_log(-x, y);
  333         // ^This is not technically correct, but it allows
  334         // functions such as cbrt(x^5), that is, x^(5/3),
  335         // to be evaluated when x is negative.
  336         // It is too complicated (and slow) to test whether y
  337         // is a formed from a ratio of an integer to an odd integer.
  338         // (And due to floating point inaccuracy, pointless too.)
  339         // For example, x^1.30769230769... is
  340         // actually x^(17/13), i.e. (x^17) ^ (1/13).
  341         // (-5)^(17/13) gives us now -8.204227562330453.
  342         // To see whether the result is right, we can test the given
  343         // root: (-8.204227562330453)^13 gives us the value of (-5)^17,
  344         // which proves that the expression was correct.
  345         //
  346         // The y*16 check prevents e.g. (-4)^(3/2) from being calculated,
  347         // as it would confuse functioninfo when pow() returns no error
  348         // but sqrt() does when the formula is converted into sqrt(x)*x.
  349         //
  350         // The errors in this approach are:
  351         //     (-2)^sqrt(2) should produce NaN
  352         //                  or actually sqrt(2)I + 2^sqrt(2),
  353         //                  produces -(2^sqrt(2)) instead.
  354         //                  (Impact: Neglible)
  355         // Thus, at worst, we're changing a NaN (or complex)
  356         // result into a negative real number result.
  357     }
  358     else
  359     {
  360         // y is negative. Utilize the x^y = 1/(x^-y) identity.
  361         if(x > ValueT(0)) return fp_pow_with_exp_log(ValueT(1) / x, -y);
  362         if(x < ValueT(0))
  363         {
  364             if(!isInteger(y*ValueT(-16)))
  365                 return -fp_pow_with_exp_log(ValueT(-1) / x, -y);
  366             // ^ See comment above.
  367         }
  368         // Remaining case: 0.0 ^ negative number
  369     }
  370     // This is reached when:
  371     //      x=0, and y<0
  372     //      x<0, and y*16 is either positive or negative integer
  373     // It is used for producing error values and as a safe fallback.
  374     return fp_pow_base(x, y);
  375 }
  376 #endif
  377 
  378 
  379 //=========================================================================
  380 // Elementary (atom) parsing functions
  381 //=========================================================================
  382 namespace
  383 {
  384     const unsigned FP_ParamGuardMask = 1U << (sizeof(unsigned) * 8u - 1u);
  385     // ^ This mask is used to prevent cFetch/other opcode's parameters
  386     //   from being confused into opcodes or variable indices within the
  387     //   bytecode optimizer. Because the way it is tested in bytecoderules.dat
  388     //   for speed reasons, it must also be the sign-bit of the "int" datatype.
  389     //   Perhaps an "assert(IsVarOpcode(X | FP_ParamGuardMask) == false)"
  390     //   might be justified to put somewhere in the code, just in case?
  391 
  392 
  393     /* Reads an UTF8-encoded sequence which forms a valid identifier name from
  394        the given input string and returns its length. If bit 31 is set, the
  395        return value also contains the internal function opcode (defined in
  396        fptypes.hh) that matches the name.
  397     */
  398     unsigned readIdentifierCommon(const char* input)
  399     {
  400         /* Assuming unsigned = 32 bits:
  401               76543210 76543210 76543210 76543210
  402            Return value if built-in function:
  403               1PPPPPPP PPPPPPPP LLLLLLLL LLLLLLLL
  404                 P = function opcode      (15 bits)
  405                 L = function name length (16 bits)
  406            Return value if not built-in function:
  407               0LLLLLLL LLLLLLLL LLLLLLLL LLLLLLLL
  408                 L = identifier length (31 bits)
  409            If unsigned has more than 32 bits, the other
  410            higher order bits are to be assumed zero.
  411         */
  412 #include "extrasrc/fp_identifier_parser.inc"
  413         return 0;
  414     }
  415 
  416     template<typename Value_t>
  417     inline unsigned readIdentifier(const char* input)
  418     {
  419         const unsigned value = readIdentifierCommon(input);
  420         if( (value & 0x80000000U) != 0) // Function?
  421         {
  422             // Verify that the function actually exists for this datatype
  423             if(IsIntType<Value_t>::result
  424             && !Functions[(value >> 16) & 0x7FFF].okForInt())
  425             {
  426                 // If it does not exist, return it as an identifier instead
  427                 return value & 0xFFFFu;
  428             }
  429             if(!IsComplexType<Value_t>::result
  430             && Functions[(value >> 16) & 0x7FFF].complexOnly())
  431             {
  432                 // If it does not exist, return it as an identifier instead
  433                 return value & 0xFFFFu;
  434             }
  435         }
  436         return value;
  437     }
  438 
  439     // Returns true if the entire string is a valid identifier
  440     template<typename Value_t>
  441     bool containsOnlyValidIdentifierChars(const std::string& name)
  442     {
  443         if(name.empty()) return false;
  444         return readIdentifier<Value_t>(name.c_str()) == unsigned(name.size());
  445     }
  446 
  447 
  448     // -----------------------------------------------------------------------
  449     // Wrappers for strto... functions
  450     // -----------------------------------------------------------------------
  451     template<typename Value_t>
  452     inline Value_t fp_parseLiteral(const char* str, char** endptr)
  453     {
  454         return std::strtod(str, endptr);
  455     }
  456 
  457 #if defined(FP_USE_STRTOLD) || defined(FP_SUPPORT_CPLUSPLUS11_MATH_FUNCS)
  458     template<>
  459     inline long double fp_parseLiteral<long double>(const char* str,
  460                                                     char** endptr)
  461     {
  462         using namespace std; // Just in case strtold() is not inside std::
  463         return strtold(str, endptr);
  464     }
  465 #endif
  466 
  467 #ifdef FP_SUPPORT_LONG_INT_TYPE
  468     template<>
  469     inline long fp_parseLiteral<long>(const char* str, char** endptr)
  470     {
  471         return std::strtol(str, endptr, 10);
  472     }
  473 #endif
  474 
  475 #ifdef FP_SUPPORT_COMPLEX_NUMBERS
  476     template<typename T>
  477     inline std::complex<T> fp_parseComplexLiteral(const char* str,
  478                                                   char** endptr)
  479     {
  480         T result = fp_parseLiteral<T> (str,endptr);
  481         const char* end = *endptr;
  482         if( (*end == 'i'  || *end == 'I')
  483         &&  !std::isalnum(end[1]) )
  484         {
  485             ++*endptr;
  486             return std::complex<T> (T(), result);
  487         }
  488         return std::complex<T> (result, T());
  489     }
  490 #endif
  491 
  492 #ifdef FP_SUPPORT_COMPLEX_DOUBLE_TYPE
  493     template<>
  494     inline std::complex<double> fp_parseLiteral<std::complex<double> >
  495     (const char* str, char** endptr)
  496     {
  497         return fp_parseComplexLiteral<double> (str,endptr);
  498     }
  499 #endif
  500 
  501 #ifdef FP_SUPPORT_COMPLEX_FLOAT_TYPE
  502     template<>
  503     inline std::complex<float> fp_parseLiteral<std::complex<float> >
  504     (const char* str, char** endptr)
  505     {
  506         return fp_parseComplexLiteral<float> (str,endptr);
  507     }
  508 #endif
  509 
  510 #ifdef FP_SUPPORT_COMPLEX_LONG_DOUBLE_TYPE
  511     template<>
  512     inline std::complex<long double> fp_parseLiteral<std::complex<long double> >
  513     (const char* str, char** endptr)
  514     {
  515         return fp_parseComplexLiteral<long double> (str,endptr);
  516     }
  517 #endif
  518 
  519     // -----------------------------------------------------------------------
  520     // Hexadecimal floating point literal parsing
  521     // -----------------------------------------------------------------------
  522     inline int testXdigit(unsigned c)
  523     {
  524         if((c-'0') < 10u) return c&15; // 0..9
  525         if(((c|0x20)-'a') < 6u) return 9+(c&15); // A..F or a..f
  526         return -1; // Not a hex digit
  527     }
  528 
  529     template<typename elem_t, unsigned n_limbs, unsigned limb_bits>
  530     inline void addXdigit(elem_t* buffer, unsigned nibble)
  531     {
  532         for(unsigned p=0; p<n_limbs; ++p)
  533         {
  534             unsigned carry = unsigned( buffer[p] >> elem_t(limb_bits-4) );
  535             buffer[p] = (buffer[p] << 4) | nibble;
  536             nibble = carry;
  537         }
  538     }
  539 
  540     template<typename Value_t>
  541     Value_t parseHexLiteral(const char* str, char** endptr)
  542     {
  543         const unsigned bits_per_char = 8;
  544 
  545         const int MantissaBits =
  546             std::numeric_limits<Value_t>::radix == 2
  547             ? std::numeric_limits<Value_t>::digits
  548             : (((sizeof(Value_t) * bits_per_char) &~ 3) - 4);
  549 
  550         typedef unsigned long elem_t;
  551         const int ExtraMantissaBits = 4 + ((MantissaBits+3)&~3); // Store one digit more for correct rounding
  552         const unsigned limb_bits = sizeof(elem_t) * bits_per_char;
  553         const unsigned n_limbs   = (ExtraMantissaBits + limb_bits-1) / limb_bits;
  554         elem_t mantissa_buffer[n_limbs] = { 0 };
  555 
  556         int n_mantissa_bits = 0; // Track the number of bits
  557         int exponent = 0; // The exponent that will be used to multiply the mantissa
  558         // Read integer portion
  559         while(true)
  560         {
  561             int xdigit = testXdigit(*str);
  562             if(xdigit < 0) break;
  563             addXdigit<elem_t,n_limbs,limb_bits> (mantissa_buffer, xdigit);
  564             ++str;
  565 
  566             n_mantissa_bits += 4;
  567             if(n_mantissa_bits >= ExtraMantissaBits)
  568             {
  569                 // Exhausted the precision. Parse the rest (until exponent)
  570                 // normally but ignore the actual digits.
  571                 for(; testXdigit(*str) >= 0; ++str)
  572                     exponent += 4;
  573                 // Read but ignore decimals
  574                 if(*str == '.')
  575                     for(++str; testXdigit(*str) >= 0; ++str)
  576                         {}
  577                 goto read_exponent;
  578             }
  579         }
  580         // Read decimals
  581         if(*str == '.')
  582             for(++str; ; )
  583             {
  584                 int xdigit = testXdigit(*str);
  585                 if(xdigit < 0) break;
  586                 addXdigit<elem_t,n_limbs,limb_bits> (mantissa_buffer, xdigit);
  587                 ++str;
  588 
  589                 exponent -= 4;
  590                 n_mantissa_bits += 4;
  591                 if(n_mantissa_bits >= ExtraMantissaBits)
  592                 {
  593                     // Exhausted the precision. Skip the rest
  594                     // of the decimals, until the exponent.
  595                     while(testXdigit(*str) >= 0)
  596                         ++str;
  597                     break;
  598                 }
  599             }
  600 
  601         // Read exponent
  602     read_exponent:
  603         if(*str == 'p' || *str == 'P')
  604         {
  605             const char* str2 = str+1;
  606             long p_exponent = strtol(str2, const_cast<char**> (&str2), 10);
  607             if(str2 != str+1 && p_exponent == (long)(int)p_exponent)
  608             {
  609                 exponent += (int)p_exponent;
  610                 str = str2;
  611             }
  612         }
  613 
  614         if(endptr) *endptr = const_cast<char*> (str);
  615 
  616         Value_t result = std::ldexp(Value_t(mantissa_buffer[0]), exponent);
  617         for(unsigned p=1; p<n_limbs; ++p)
  618         {
  619             exponent += limb_bits;
  620             result += ldexp(Value_t(mantissa_buffer[p]), exponent);
  621         }
  622         return result;
  623     }
  624 
  625 #ifdef FP_SUPPORT_LONG_INT_TYPE
  626     template<>
  627     long parseHexLiteral<long>(const char* str, char** endptr)
  628     {
  629         return std::strtol(str, endptr, 16);
  630     }
  631 #endif
  632 
  633 #ifdef FP_SUPPORT_COMPLEX_DOUBLE_TYPE
  634     template<>
  635     std::complex<double>
  636     parseHexLiteral<std::complex<double> >(const char* str, char** endptr)
  637     {
  638         return parseHexLiteral<double> (str, endptr);
  639     }
  640 #endif
  641 
  642 #ifdef FP_SUPPORT_COMPLEX_FLOAT_TYPE
  643     template<>
  644     std::complex<float>
  645     parseHexLiteral<std::complex<float> >(const char* str, char** endptr)
  646     {
  647         return parseHexLiteral<float> (str, endptr);
  648     }
  649 #endif
  650 
  651 #ifdef FP_SUPPORT_COMPLEX_LONG_DOUBLE_TYPE
  652     template<>
  653     std::complex<long double>
  654     parseHexLiteral<std::complex<long double> >(const char* str, char** endptr)
  655     {
  656         return parseHexLiteral<long double> (str, endptr);
  657     }
  658 #endif
  659 }
  660 
  661 //=========================================================================
  662 // Utility functions
  663 //=========================================================================
  664 namespace
  665 {
  666     // -----------------------------------------------------------------------
  667     // Add a new identifier to the specified identifier map
  668     // -----------------------------------------------------------------------
  669     // Return value will be false if the name already existed
  670     template<typename Value_t>
  671     bool addNewNameData(NamePtrsMap<Value_t>& namePtrs,
  672                         std::pair<NamePtr, NameData<Value_t> >& newName,
  673                         bool isVar)
  674     {
  675         typename NamePtrsMap<Value_t>::iterator nameIter =
  676             namePtrs.lower_bound(newName.first);
  677 
  678         if(nameIter != namePtrs.end() && newName.first == nameIter->first)
  679         {
  680             // redefining a var is not allowed.
  681             if(isVar) return false;
  682 
  683             // redefining other tokens is allowed, if the type stays the same.
  684             if(nameIter->second.type != newName.second.type)
  685                 return false;
  686 
  687             // update the data
  688             nameIter->second = newName.second;
  689             return true;
  690         }
  691 
  692         if(!isVar)
  693         {
  694             // Allocate a copy of the name (pointer stored in the map key)
  695             // However, for VARIABLEs, the pointer points to VariableString,
  696             // which is managed separately. Thusly, only done when !IsVar.
  697             char* namebuf = new char[newName.first.nameLength];
  698             memcpy(namebuf, newName.first.name, newName.first.nameLength);
  699             newName.first.name = namebuf;
  700         }
  701 
  702         namePtrs.insert(nameIter, newName);
  703         return true;
  704     }
  705 }
  706 
  707 
  708 //=========================================================================
  709 // Data struct implementation
  710 //=========================================================================
  711 template<typename Value_t>
  712 FunctionParserBase<Value_t>::Data::Data():
  713     mReferenceCounter(1),
  714     mDelimiterChar(0),
  715     mParseErrorType(NO_FUNCTION_PARSED_YET),
  716     mEvalErrorType(0),
  717     mUseDegreeConversion(false),
  718     mErrorLocation(0),
  719     mVariablesAmount(0),
  720     mStackSize(0)
  721 {}
  722 
  723 template<typename Value_t>
  724 FunctionParserBase<Value_t>::Data::Data(const Data& rhs):
  725     mReferenceCounter(0),
  726     mDelimiterChar(rhs.mDelimiterChar),
  727     mParseErrorType(rhs.mParseErrorType),
  728     mEvalErrorType(rhs.mEvalErrorType),
  729     mUseDegreeConversion(rhs.mUseDegreeConversion),
  730     mErrorLocation(rhs.mErrorLocation),
  731     mVariablesAmount(rhs.mVariablesAmount),
  732     mVariablesString(rhs.mVariablesString),
  733     mNamePtrs(),
  734     mFuncPtrs(rhs.mFuncPtrs),
  735     mFuncParsers(rhs.mFuncParsers),
  736     mByteCode(rhs.mByteCode),
  737     mImmed(rhs.mImmed),
  738 #ifndef FP_USE_THREAD_SAFE_EVAL
  739     mStack(rhs.mStackSize),
  740 #endif
  741     mStackSize(rhs.mStackSize)
  742 {
  743     for(typename NamePtrsMap<Value_t>::const_iterator i = rhs.mNamePtrs.begin();
  744         i != rhs.mNamePtrs.end();
  745         ++i)
  746     {
  747         if(i->second.type == NameData<Value_t>::VARIABLE)
  748         {
  749             const std::size_t variableStringOffset =
  750                 i->first.name - rhs.mVariablesString.c_str();
  751             std::pair<NamePtr, NameData<Value_t> > tmp
  752                 (NamePtr(&mVariablesString[variableStringOffset],
  753                          i->first.nameLength),
  754                  i->second);
  755             mNamePtrs.insert(mNamePtrs.end(), tmp);
  756         }
  757         else
  758         {
  759             std::pair<NamePtr, NameData<Value_t> > tmp
  760                 (NamePtr(new char[i->first.nameLength], i->first.nameLength),
  761                  i->second );
  762             memcpy(const_cast<char*>(tmp.first.name), i->first.name,
  763                    tmp.first.nameLength);
  764             mNamePtrs.insert(mNamePtrs.end(), tmp);
  765         }
  766     }
  767 }
  768 
  769 template<typename Value_t>
  770 FunctionParserBase<Value_t>::Data::~Data()
  771 {
  772     for(typename NamePtrsMap<Value_t>::iterator i = mNamePtrs.begin();
  773         i != mNamePtrs.end();
  774         ++i)
  775     {
  776         if(i->second.type != NameData<Value_t>::VARIABLE)
  777             delete[] i->first.name;
  778     }
  779 }
  780 
  781 template<typename Value_t>
  782 void FunctionParserBase<Value_t>::incFuncWrapperRefCount
  783 (FunctionWrapper* wrapper)
  784 {
  785     ++wrapper->mReferenceCount;
  786 }
  787 
  788 template<typename Value_t>
  789 unsigned FunctionParserBase<Value_t>::decFuncWrapperRefCount
  790 (FunctionWrapper* wrapper)
  791 {
  792     return --wrapper->mReferenceCount;
  793 }
  794 
  795 template<typename Value_t>
  796 FunctionParserBase<Value_t>::Data::FuncWrapperPtrData::FuncWrapperPtrData():
  797     mRawFuncPtr(0), mFuncWrapperPtr(0), mParams(0)
  798 {}
  799 
  800 template<typename Value_t>
  801 FunctionParserBase<Value_t>::Data::FuncWrapperPtrData::~FuncWrapperPtrData()
  802 {
  803     if(mFuncWrapperPtr &&
  804        FunctionParserBase::decFuncWrapperRefCount(mFuncWrapperPtr) == 0)
  805         delete mFuncWrapperPtr;
  806 }
  807 
  808 template<typename Value_t>
  809 FunctionParserBase<Value_t>::Data::FuncWrapperPtrData::FuncWrapperPtrData
  810 (const FuncWrapperPtrData& rhs):
  811     mRawFuncPtr(rhs.mRawFuncPtr),
  812     mFuncWrapperPtr(rhs.mFuncWrapperPtr),
  813     mParams(rhs.mParams)
  814 {
  815     if(mFuncWrapperPtr)
  816         FunctionParserBase::incFuncWrapperRefCount(mFuncWrapperPtr);
  817 }
  818 
  819 template<typename Value_t>
  820 typename FunctionParserBase<Value_t>::Data::FuncWrapperPtrData&
  821 FunctionParserBase<Value_t>::Data::FuncWrapperPtrData::operator=
  822 (const FuncWrapperPtrData& rhs)
  823 {
  824     if(&rhs != this)
  825     {
  826         if(mFuncWrapperPtr &&
  827            FunctionParserBase::decFuncWrapperRefCount(mFuncWrapperPtr) == 0)
  828             delete mFuncWrapperPtr;
  829         mRawFuncPtr = rhs.mRawFuncPtr;
  830         mFuncWrapperPtr = rhs.mFuncWrapperPtr;
  831         mParams = rhs.mParams;
  832         if(mFuncWrapperPtr)
  833             FunctionParserBase::incFuncWrapperRefCount(mFuncWrapperPtr);
  834     }
  835     return *this;
  836 }
  837 
  838 
  839 //=========================================================================
  840 // FunctionParser constructors, destructor and assignment
  841 //=========================================================================
  842 template<typename Value_t>
  843 FunctionParserBase<Value_t>::FunctionParserBase():
  844     mData(new Data),
  845     mStackPtr(0)
  846 {
  847 }
  848 
  849 template<typename Value_t>
  850 FunctionParserBase<Value_t>::~FunctionParserBase()
  851 {
  852     if(--(mData->mReferenceCounter) == 0)
  853         delete mData;
  854 }
  855 
  856 template<typename Value_t>
  857 FunctionParserBase<Value_t>::FunctionParserBase(const FunctionParserBase& cpy):
  858     mData(cpy.mData),
  859     mStackPtr(0)
  860 {
  861     ++(mData->mReferenceCounter);
  862 }
  863 
  864 template<typename Value_t>
  865 FunctionParserBase<Value_t>&
  866 FunctionParserBase<Value_t>::operator=(const FunctionParserBase& cpy)
  867 {
  868     if(mData != cpy.mData)
  869     {
  870         if(--(mData->mReferenceCounter) == 0) delete mData;
  871 
  872         mData = cpy.mData;
  873         ++(mData->mReferenceCounter);
  874     }
  875     return *this;
  876 }
  877 
  878 template<typename Value_t>
  879 typename FunctionParserBase<Value_t>::Data*
  880 FunctionParserBase<Value_t>::getParserData()
  881 {
  882     return mData;
  883 }
  884 
  885 template<typename Value_t>
  886 void FunctionParserBase<Value_t>::setDelimiterChar(char c)
  887 {
  888     mData->mDelimiterChar = c;
  889 }
  890 
  891 
  892 //---------------------------------------------------------------------------
  893 // Copy-on-write method
  894 //---------------------------------------------------------------------------
  895 template<typename Value_t>
  896 void FunctionParserBase<Value_t>::CopyOnWrite()
  897 {
  898     if(mData->mReferenceCounter > 1)
  899     {
  900         Data* oldData = mData;
  901         mData = new Data(*oldData);
  902         --(oldData->mReferenceCounter);
  903         mData->mReferenceCounter = 1;
  904     }
  905 }
  906 
  907 template<typename Value_t>
  908 void FunctionParserBase<Value_t>::CopyOnWrite2()
  909 {
  910     if(mData->mReferenceCounter > 1)
  911     {
  912         Data* oldData = mData;
  913         mData = new Data(*oldData);
  914         //--(oldData->mReferenceCounter);
  915         //mData->mReferenceCounter = 1;
  916         for (unsigned i=0; i<oldData->mFuncParsers.size(); i++)
  917         {
  918             mData->mFuncParsers[i].mParserPtr = new FunctionParserBase(*(oldData->mFuncParsers[i].mParserPtr));
  919             mData->mFuncParsers[i].mParams = oldData->mFuncParsers[i].mParams;
  920         }
  921         mData->mByteCode = oldData->mByteCode;
  922 
  923     }
  924 }
  925 
  926 template<typename Value_t>
  927 void FunctionParserBase<Value_t>::ForceDeepCopy()
  928 {
  929     CopyOnWrite2();
  930 }
  931 
  932 //=========================================================================
  933 // Epsilon
  934 //=========================================================================
  935 template<typename Value_t>
  936 Value_t FunctionParserBase<Value_t>::epsilon()
  937 {
  938     return Epsilon<Value_t>::value;
  939 }
  940 
  941 template<typename Value_t>
  942 void FunctionParserBase<Value_t>::setEpsilon(Value_t value)
  943 {
  944     Epsilon<Value_t>::value = value;
  945 }
  946 
  947 //=========================================================================
  948 // User-defined identifier addition functions
  949 //=========================================================================
  950 template<typename Value_t>
  951 bool FunctionParserBase<Value_t>::AddConstant(const std::string& name, Value_t value)
  952 {
  953     if(!containsOnlyValidIdentifierChars<Value_t>(name)) return false;
  954 
  955     CopyOnWrite();
  956     std::pair<NamePtr, NameData<Value_t> > newName
  957         (NamePtr(name.data(), unsigned(name.size())),
  958          NameData<Value_t>(NameData<Value_t>::CONSTANT, value));
  959 
  960     return addNewNameData(mData->mNamePtrs, newName, false);
  961 }
  962 
  963 template<typename Value_t>
  964 bool FunctionParserBase<Value_t>::AddUnit(const std::string& name,
  965                                           Value_t value)
  966 {
  967     if(!containsOnlyValidIdentifierChars<Value_t>(name)) return false;
  968 
  969     CopyOnWrite();
  970     std::pair<NamePtr, NameData<Value_t> > newName
  971         (NamePtr(name.data(), unsigned(name.size())),
  972          NameData<Value_t>(NameData<Value_t>::UNIT, value));
  973     return addNewNameData(mData->mNamePtrs, newName, false);
  974 }
  975 
  976 template<typename Value_t>
  977 bool FunctionParserBase<Value_t>::AddFunction
  978 (const std::string& name, FunctionPtr ptr, unsigned paramsAmount)
  979 {
  980     if(!containsOnlyValidIdentifierChars<Value_t>(name)) return false;
  981 
  982     CopyOnWrite();
  983     std::pair<NamePtr, NameData<Value_t> > newName
  984         (NamePtr(name.data(), unsigned(name.size())),
  985          NameData<Value_t>(NameData<Value_t>::FUNC_PTR,
  986                            unsigned(mData->mFuncPtrs.size())));
  987 
  988     const bool success = addNewNameData(mData->mNamePtrs, newName, false);
  989     if(success)
  990     {
  991         mData->mFuncPtrs.push_back(typename Data::FuncWrapperPtrData());
  992         mData->mFuncPtrs.back().mRawFuncPtr = ptr;
  993         mData->mFuncPtrs.back().mParams = paramsAmount;
  994     }
  995     return success;
  996 }
  997 
  998 template<typename Value_t>
  999 bool FunctionParserBase<Value_t>::addFunctionWrapperPtr
 1000 (const std::string& name, FunctionWrapper* wrapper, unsigned paramsAmount)
 1001 {
 1002     if(!AddFunction(name, FunctionPtr(0), paramsAmount)) return false;
 1003     mData->mFuncPtrs.back().mFuncWrapperPtr = wrapper;
 1004     return true;
 1005 }
 1006 
 1007 template<typename Value_t>
 1008 typename FunctionParserBase<Value_t>::FunctionWrapper*
 1009 FunctionParserBase<Value_t>::GetFunctionWrapper(const std::string& name)
 1010 {
 1011     CopyOnWrite();
 1012     NamePtr namePtr(name.data(), unsigned(name.size()));
 1013 
 1014     typename NamePtrsMap<Value_t>::iterator nameIter =
 1015         mData->mNamePtrs.find(namePtr);
 1016 
 1017     if(nameIter != mData->mNamePtrs.end() &&
 1018        nameIter->second.type == NameData<Value_t>::FUNC_PTR)
 1019     {
 1020         return mData->mFuncPtrs[nameIter->second.index].mFuncWrapperPtr;
 1021     }
 1022     return 0;
 1023 }
 1024 
 1025 template<typename Value_t>
 1026 bool FunctionParserBase<Value_t>::CheckRecursiveLinking
 1027 (const FunctionParserBase* fp) const
 1028 {
 1029     if(fp == this) return true;
 1030     for(unsigned i = 0; i < fp->mData->mFuncParsers.size(); ++i)
 1031         if(CheckRecursiveLinking(fp->mData->mFuncParsers[i].mParserPtr))
 1032             return true;
 1033     return false;
 1034 }
 1035 
 1036 template<typename Value_t>
 1037 bool FunctionParserBase<Value_t>::AddFunction(const std::string& name,
 1038                                               FunctionParserBase& fp)
 1039 {
 1040     if(!containsOnlyValidIdentifierChars<Value_t>(name) ||
 1041        CheckRecursiveLinking(&fp))
 1042         return false;
 1043 
 1044     CopyOnWrite();
 1045     std::pair<NamePtr, NameData<Value_t> > newName
 1046         (NamePtr(name.data(), unsigned(name.size())),
 1047          NameData<Value_t>(NameData<Value_t>::PARSER_PTR,
 1048                            unsigned(mData->mFuncParsers.size())));
 1049 
 1050     const bool success = addNewNameData(mData->mNamePtrs, newName, false);
 1051     if(success)
 1052     {
 1053         mData->mFuncParsers.push_back(typename Data::FuncParserPtrData());
 1054         mData->mFuncParsers.back().mParserPtr = &fp;
 1055         mData->mFuncParsers.back().mParams = fp.mData->mVariablesAmount;
 1056     }
 1057     return success;
 1058 }
 1059 
 1060 template<typename Value_t>
 1061 bool FunctionParserBase<Value_t>::RemoveIdentifier(const std::string& name)
 1062 {
 1063     CopyOnWrite();
 1064 
 1065     NamePtr namePtr(name.data(), unsigned(name.size()));
 1066 
 1067     typename NamePtrsMap<Value_t>::iterator nameIter =
 1068         mData->mNamePtrs.find(namePtr);
 1069 
 1070     if(nameIter != mData->mNamePtrs.end())
 1071     {
 1072         if(nameIter->second.type == NameData<Value_t>::VARIABLE)
 1073         {
 1074             // Illegal attempt to delete variables
 1075             return false;
 1076         }
 1077         delete[] nameIter->first.name;
 1078         mData->mNamePtrs.erase(nameIter);
 1079         return true;
 1080     }
 1081     return false;
 1082 }
 1083 
 1084 //=========================================================================
 1085 // Function parsing
 1086 //=========================================================================
 1087 namespace
 1088 {
 1089     // Error messages returned by ErrorMsg():
 1090     const char* const ParseErrorMessage[]=
 1091     {
 1092         "Syntax error",                             // 0
 1093         "Mismatched parenthesis",                   // 1
 1094         "Missing ')'",                              // 2
 1095         "Empty parentheses",                        // 3
 1096         "Syntax error: Operator expected",          // 4
 1097         "Not enough memory",                        // 5
 1098         "An unexpected error occurred. Please make a full bug report "
 1099         "to the author",                            // 6
 1100         "Syntax error in parameter 'Vars' given to "
 1101         "FunctionParser::Parse()",                  // 7
 1102         "Illegal number of parameters to function", // 8
 1103         "Syntax error: Premature end of string",    // 9
 1104         "Syntax error: Expecting ( after function", // 10
 1105         "Syntax error: Unknown identifier",         // 11
 1106         "(No function has been parsed yet)",
 1107         ""
 1108     };
 1109 
 1110     template<typename Value_t>
 1111     inline typename FunctionParserBase<Value_t>::ParseErrorType
 1112     noCommaError(char c)
 1113     {
 1114         return c == ')' ?
 1115             FunctionParserBase<Value_t>::ILL_PARAMS_AMOUNT :
 1116             FunctionParserBase<Value_t>::SYNTAX_ERROR;
 1117     }
 1118 
 1119     template<typename Value_t>
 1120     inline typename FunctionParserBase<Value_t>::ParseErrorType
 1121     noParenthError(char c)
 1122     {
 1123         return c == ',' ?
 1124             FunctionParserBase<Value_t>::ILL_PARAMS_AMOUNT :
 1125             FunctionParserBase<Value_t>::MISSING_PARENTH;
 1126     }
 1127 
 1128     template<unsigned offset>
 1129     struct IntLiteralMask
 1130     {
 1131         enum { mask =
 1132         //    (    1UL << ('-'-offset)) |
 1133             (0x3FFUL << ('0'-offset)) }; /* 0x3FF = 10 bits worth "1" */
 1134         // Note: If you change fparser to support negative numbers parsing
 1135         //       (as opposed to parsing them as cNeg followed by literal),
 1136         //       enable the '-' line above, and change the offset value
 1137         //       in BeginsLiteral() to '-' instead of '.'.
 1138     };
 1139 
 1140     template<typename Value_t, unsigned offset>
 1141     struct LiteralMask
 1142     {
 1143         enum { mask =
 1144             (    1UL << ('.'-offset)) |
 1145             IntLiteralMask<offset>::mask };
 1146     };
 1147 #ifdef FP_SUPPORT_LONG_INT_TYPE
 1148     template<unsigned offset>
 1149     struct LiteralMask<long, offset>: public IntLiteralMask<offset>
 1150     {
 1151     };
 1152 #endif
 1153 #ifdef FP_SUPPORT_GMP_INT_TYPE
 1154     template<unsigned offset>
 1155     struct LiteralMask<GmpInt, offset>: public IntLiteralMask<offset>
 1156     {
 1157     };
 1158 #endif
 1159 
 1160     template<unsigned offset>
 1161     struct SimpleSpaceMask
 1162     {
 1163         enum { mask =
 1164             (1UL << ('\r'-offset)) |
 1165             (1UL << ('\n'-offset)) |
 1166             (1UL << ('\v'-offset)) |
 1167             (1UL << ('\t'-offset)) |
 1168             (1UL << (' ' -offset)) };
 1169     };
 1170 
 1171     template<typename Value_t>
 1172     inline bool BeginsLiteral(unsigned byte)
 1173     {
 1174         enum { n = sizeof(unsigned long)>=8 ? 0 : '.' };
 1175         byte -= n;
 1176         if(byte > (unsigned char)('9'-n)) return false;
 1177         unsigned long shifted = 1UL << byte;
 1178         const unsigned long mask = LiteralMask<Value_t, n>::mask;
 1179         return (mask & shifted) != 0;
 1180     }
 1181 
 1182     template<typename CharPtr>
 1183     inline void SkipSpace(CharPtr& function)
 1184     {
 1185 /*
 1186         Space characters in unicode:
 1187 U+0020  SPACE                      Depends on font, often adjusted (see below)
 1188 U+00A0  NO-BREAK SPACE             As a space, but often not adjusted
 1189 U+2000  EN QUAD                    1 en (= 1/2 em)
 1190 U+2001  EM QUAD                    1 em (nominally, the height of the font)
 1191 U+2002  EN SPACE                   1 en (= 1/2 em)
 1192 U+2003  EM SPACE                   1 em
 1193 U+2004  THREE-PER-EM SPACE         1/3 em
 1194 U+2005  FOUR-PER-EM SPACE          1/4 em
 1195 U+2006  SIX-PER-EM SPACE           1/6 em
 1196 U+2007  FIGURE SPACE               Tabular width, the width of digits
 1197 U+2008  PUNCTUATION SPACE          The width of a period .
 1198 U+2009  THIN SPACE                 1/5 em (or sometimes 1/6 em)
 1199 U+200A  HAIR SPACE                 Narrower than THIN SPACE
 1200 U+200B  ZERO WIDTH SPACE           Nominally no width, but may expand
 1201 U+202F  NARROW NO-BREAK SPACE      Narrower than NO-BREAK SPACE (or SPACE)
 1202 U+205F  MEDIUM MATHEMATICAL SPACE  4/18 em
 1203 U+3000  IDEOGRAPHIC SPACE          The width of ideographic (CJK) characters.
 1204         Also:
 1205 U+000A  \n
 1206 U+000D  \r
 1207 U+0009  \t
 1208 U+000B  \v
 1209         As UTF-8 sequences:
 1210             09
 1211             0A
 1212             0B
 1213             0D
 1214             20
 1215             C2 A0
 1216             E2 80 80-8B
 1217             E2 80 AF
 1218             E2 81 9F
 1219             E3 80 80
 1220 */
 1221         while(true)
 1222         {
 1223             enum { n = sizeof(unsigned long)>=8 ? 0 : '\t' };
 1224             typedef signed char schar;
 1225             unsigned byte = (unsigned char)*function;
 1226             byte -= n;
 1227             // ^Note: values smaller than n intentionally become
 1228             //        big values here due to integer wrap. The
 1229             //        comparison below thus excludes them, making
 1230             //        the effective range 0x09..0x20 (32-bit)
 1231             //        or 0x00..0x20 (64-bit) within the if-clause.
 1232             if(byte <= (unsigned char)(' '-n))
 1233             {
 1234                 unsigned long shifted = 1UL << byte;
 1235                 const unsigned long mask = SimpleSpaceMask<n>::mask;
 1236                 if(mask & shifted)
 1237                     { ++function; continue; } // \r, \n, \t, \v and space
 1238                 break;
 1239             }
 1240             if(likely(byte < 0xC2-n)) break;
 1241 
 1242             if(byte == 0xC2-n && function[1] == char(0xA0))
 1243                 { function += 2; continue; } // U+00A0
 1244             if(byte == 0xE3-n &&
 1245                function[1] == char(0x80) && function[2] == char(0x80))
 1246                 { function += 3; continue; } // U+3000
 1247             if(byte == 0xE2-n)
 1248             {
 1249                 if(function[1] == char(0x81))
 1250                 {
 1251                     if(function[2] != char(0x9F)) break;
 1252                     function += 3; // U+205F
 1253                     continue;
 1254                 }
 1255                 if(function[1] == char(0x80))
 1256                 if(function[2] == char(0xAF) || // U+202F
 1257                    schar(function[2]) <= schar(0x8B) // U+2000..U+200B
 1258                   )
 1259                 {
 1260                     function += 3;
 1261                     continue;
 1262                 }
 1263             }
 1264             break;
 1265         } // while(true)
 1266     } // SkipSpace(CharPtr& function)
 1267 }
 1268 
 1269 // ---------------------------------------------------------------------------
 1270 // Return parse error message
 1271 // ---------------------------------------------------------------------------
 1272 template<typename Value_t>
 1273 const char* FunctionParserBase<Value_t>::ErrorMsg() const
 1274 {
 1275     return ParseErrorMessage[mData->mParseErrorType];
 1276 }
 1277 
 1278 template<typename Value_t>
 1279 typename FunctionParserBase<Value_t>::ParseErrorType
 1280 FunctionParserBase<Value_t>::GetParseErrorType() const
 1281 {
 1282     return mData->mParseErrorType;
 1283 }
 1284 
 1285 template<typename Value_t>
 1286 int FunctionParserBase<Value_t>::EvalError() const
 1287 {
 1288     return mData->mEvalErrorType;
 1289 }
 1290 
 1291 // ---------------------------------------------------------------------------
 1292 // Parse variables
 1293 // ---------------------------------------------------------------------------
 1294 template<typename Value_t>
 1295 bool FunctionParserBase<Value_t>::ParseVariables
 1296 (const std::string& inputVarString)
 1297 {
 1298     if(mData->mVariablesString == inputVarString) return true;
 1299 
 1300     /* Delete existing variables from mNamePtrs */
 1301     for(typename NamePtrsMap<Value_t>::iterator i =
 1302             mData->mNamePtrs.begin();
 1303         i != mData->mNamePtrs.end(); )
 1304     {
 1305         if(i->second.type == NameData<Value_t>::VARIABLE)
 1306         {
 1307             typename NamePtrsMap<Value_t>::iterator j (i);
 1308             ++i;
 1309             mData->mNamePtrs.erase(j);
 1310         }
 1311         else ++i;
 1312     }
 1313     mData->mVariablesString = inputVarString;
 1314     const std::string& vars = mData->mVariablesString;
 1315     const unsigned len = unsigned(vars.size());
 1316     unsigned varNumber = VarBegin;
 1317     const char* beginPtr = vars.c_str();
 1318     const char* finalPtr = beginPtr + len;
 1319     while(beginPtr < finalPtr)
 1320     {
 1321         SkipSpace(beginPtr);
 1322         unsigned nameLength = readIdentifier<Value_t>(beginPtr);
 1323         if(nameLength == 0 || (nameLength & 0x80000000U)) return false;
 1324         const char* endPtr = beginPtr + nameLength;
 1325         SkipSpace(endPtr);
 1326         if(endPtr != finalPtr && *endPtr != ',') return false;
 1327         std::pair<NamePtr, NameData<Value_t> > newName
 1328             (NamePtr(beginPtr, nameLength),
 1329              NameData<Value_t>(NameData<Value_t>::VARIABLE, varNumber++));
 1330         if(!addNewNameData(mData->mNamePtrs, newName, true))
 1331         {
 1332             return false;
 1333         }
 1334         beginPtr = endPtr + 1;
 1335     }
 1336     mData->mVariablesAmount = varNumber - VarBegin;
 1337     return true;
 1338 }
 1339 
 1340 // ---------------------------------------------------------------------------
 1341 // Parse() public interface functions
 1342 // ---------------------------------------------------------------------------
 1343 template<typename Value_t>
 1344 int FunctionParserBase<Value_t>::Parse(const char* Function,
 1345                                        const std::string& Vars,
 1346                                        bool useDegrees)
 1347 {
 1348     CopyOnWrite();
 1349 
 1350     if(!ParseVariables(Vars))
 1351     {
 1352         mData->mParseErrorType = INVALID_VARS;
 1353         return int(strlen(Function));
 1354     }
 1355 
 1356     return ParseFunction(Function, useDegrees);
 1357 }
 1358 
 1359 template<typename Value_t>
 1360 int FunctionParserBase<Value_t>::Parse(const std::string& Function,
 1361                                        const std::string& Vars,
 1362                                        bool useDegrees)
 1363 {
 1364     CopyOnWrite();
 1365 
 1366     if(!ParseVariables(Vars))
 1367     {
 1368         mData->mParseErrorType = INVALID_VARS;
 1369         return int(Function.size());
 1370     }
 1371 
 1372     return ParseFunction(Function.c_str(), useDegrees);
 1373 }
 1374 
 1375 
 1376 // ---------------------------------------------------------------------------
 1377 // Main parsing function
 1378 // ---------------------------------------------------------------------------
 1379 template<typename Value_t>
 1380 int FunctionParserBase<Value_t>::ParseFunction(const char* function,
 1381                                                bool useDegrees)
 1382 {
 1383     mData->mUseDegreeConversion = useDegrees;
 1384     mData->mParseErrorType = FP_NO_ERROR;
 1385 
 1386     mData->mInlineVarNames.clear();
 1387     mData->mByteCode.clear(); mData->mByteCode.reserve(128);
 1388     mData->mImmed.clear(); mData->mImmed.reserve(128);
 1389     mData->mStackSize = mStackPtr = 0;
 1390 
 1391     mData->mHasByteCodeFlags = false;
 1392 
 1393     const char* ptr = Compile(function);
 1394     mData->mInlineVarNames.clear();
 1395 
 1396     if(mData->mHasByteCodeFlags)
 1397     {
 1398         for(unsigned i = unsigned(mData->mByteCode.size()); i-- > 0; )
 1399             mData->mByteCode[i] &= ~FP_ParamGuardMask;
 1400     }
 1401 
 1402     if(mData->mParseErrorType != FP_NO_ERROR)
 1403         return int(mData->mErrorLocation - function);
 1404 
 1405     assert(ptr); // Should never be null at this point. It's a bug otherwise.
 1406     if(*ptr)
 1407     {
 1408         if(mData->mDelimiterChar == 0 || *ptr != mData->mDelimiterChar)
 1409             mData->mParseErrorType = EXPECT_OPERATOR;
 1410         return int(ptr - function);
 1411     }
 1412 
 1413 #ifndef FP_USE_THREAD_SAFE_EVAL
 1414     mData->mStack.resize(mData->mStackSize);
 1415     //mData->mStacki.resize(64*mData->mStackSize);
 1416 #endif
 1417 
 1418     return -1;
 1419 }
 1420 
 1421 
 1422 //=========================================================================
 1423 // Parsing and bytecode compiling functions
 1424 //=========================================================================
 1425 template<typename Value_t>
 1426 inline const char* FunctionParserBase<Value_t>::SetErrorType(ParseErrorType t,
 1427                                                              const char* pos)
 1428 {
 1429     mData->mParseErrorType = t;
 1430     mData->mErrorLocation = pos;
 1431     return 0;
 1432 }
 1433 
 1434 template<typename Value_t>
 1435 inline void FunctionParserBase<Value_t>::incStackPtr()
 1436 {
 1437     if(++mStackPtr > mData->mStackSize) ++(mData->mStackSize);
 1438 }
 1439 
 1440 namespace
 1441 {
 1442     const unsigned char powi_factor_table[128] =
 1443     {
 1444         0,1,0,0,0,0,0,0, 0, 0,0,0,0,0,0,3,/*   0 -  15 */
 1445         0,0,0,0,0,0,0,0, 0, 5,0,3,0,0,3,0,/*  16 -  31 */
 1446         0,0,0,0,0,0,0,3, 0, 0,0,0,0,5,0,0,/*  32 -  47 */
 1447         0,0,5,3,0,0,3,5, 0, 3,0,0,3,0,0,3,/*  48 -  63 */
 1448         0,0,0,0,0,0,0,0, 0, 0,0,3,0,0,3,0,/*  64 -  79 */
 1449         0,9,0,0,0,5,0,3, 0, 0,5,7,0,0,0,5,/*  80 -  95 */
 1450         0,0,0,3,5,0,3,0, 0, 3,0,0,3,0,5,3,/*  96 - 111 */
 1451         0,0,3,5,0,9,0,7, 3,11,0,3,0,5,3,0,/* 112 - 127 */
 1452     };
 1453 
 1454     inline int get_powi_factor(long abs_int_exponent)
 1455     {
 1456         if(abs_int_exponent >= int(sizeof(powi_factor_table))) return 0;
 1457         return powi_factor_table[abs_int_exponent];
 1458     }
 1459 
 1460 #if 0
 1461     int EstimatePowiComplexity(int abs_int_exponent)
 1462     {
 1463         int cost = 0;
 1464         while(abs_int_exponent > 1)
 1465         {
 1466             int factor = get_powi_factor(abs_int_exponent);
 1467             if(factor)
 1468             {
 1469                 cost += EstimatePowiComplexity(factor);
 1470                 abs_int_exponent /= factor;
 1471                 continue;
 1472             }
 1473             if(!(abs_int_exponent & 1))
 1474             {
 1475                 abs_int_exponent /= 2;
 1476                 cost += 3; // sqr
 1477             }
 1478             else
 1479             {
 1480                 cost += 4; // dup+mul
 1481                 abs_int_exponent -= 1;
 1482             }
 1483         }
 1484         return cost;
 1485     }
 1486 #endif
 1487 
 1488     bool IsEligibleIntPowiExponent(long int_exponent)
 1489     {
 1490         if(int_exponent == 0) return false;
 1491         long abs_int_exponent = int_exponent;
 1492     #if 0
 1493         int cost = 0;
 1494 
 1495         if(abs_int_exponent < 0)
 1496         {
 1497             cost += 11;
 1498             abs_int_exponent = -abs_int_exponent;
 1499         }
 1500 
 1501         cost += EstimatePowiComplexity(abs_int_exponent);
 1502 
 1503         return cost < (10*3 + 4*4);
 1504     #else
 1505         if(abs_int_exponent < 0) abs_int_exponent = -abs_int_exponent;
 1506 
 1507         return (abs_int_exponent >= 1)
 1508             && (abs_int_exponent <= 46 ||
 1509               (abs_int_exponent <= 1024 &&
 1510               (abs_int_exponent & (abs_int_exponent - 1)) == 0));
 1511     #endif
 1512     }
 1513 
 1514     /* Needed by fp_opcode_add.inc if tracing is enabled */
 1515     template<typename Value_t>
 1516     std::string findName(const NamePtrsMap<Value_t>& nameMap,
 1517                          unsigned index,
 1518                          typename NameData<Value_t>::DataType type)
 1519     {
 1520         for(typename NamePtrsMap<Value_t>::const_iterator
 1521                 iter = nameMap.begin();
 1522             iter != nameMap.end();
 1523             ++iter)
 1524         {
 1525             if(iter->second.type == type && iter->second.index == index)
 1526                 return std::string(iter->first.name,
 1527                                    iter->first.name + iter->first.nameLength);
 1528         }
 1529         return "?";
 1530     }
 1531 }
 1532 
 1533 template<typename Value_t>
 1534 inline void FunctionParserBase<Value_t>::AddImmedOpcode(Value_t value)
 1535 {
 1536     mData->mImmed.push_back(value);
 1537     mData->mByteCode.push_back(cImmed);
 1538 }
 1539 
 1540 template<typename Value_t>
 1541 inline void FunctionParserBase<Value_t>::CompilePowi(long abs_int_exponent)
 1542 {
 1543     int num_muls=0;
 1544     while(abs_int_exponent > 1)
 1545     {
 1546         long factor = get_powi_factor(abs_int_exponent);
 1547         if(factor)
 1548         {
 1549             CompilePowi(factor);
 1550             abs_int_exponent /= factor;
 1551             continue;
 1552         }
 1553         if(!(abs_int_exponent & 1))
 1554         {
 1555             abs_int_exponent /= 2;
 1556             mData->mByteCode.push_back(cSqr);
 1557             // ^ Don't put AddFunctionOpcode here,
 1558             //   it would slow down a great deal.
 1559         }
 1560         else
 1561         {
 1562             mData->mByteCode.push_back(cDup);
 1563             incStackPtr();
 1564             abs_int_exponent -= 1;
 1565             ++num_muls;
 1566         }
 1567     }
 1568     if(num_muls > 0)
 1569     {
 1570         mData->mByteCode.resize(mData->mByteCode.size()+num_muls, cMul);
 1571         mStackPtr -= num_muls;
 1572     }
 1573 }
 1574 
 1575 template<typename Value_t>
 1576 inline bool FunctionParserBase<Value_t>::TryCompilePowi(Value_t original_immed)
 1577 {
 1578     Value_t changed_immed = original_immed;
 1579     for(int sqrt_count=0; /**/; ++sqrt_count)
 1580     {
 1581         long int_exponent = makeLongInteger(changed_immed);
 1582         if(isLongInteger(changed_immed) &&
 1583            IsEligibleIntPowiExponent(int_exponent))
 1584         {
 1585             long abs_int_exponent = int_exponent;
 1586             if(abs_int_exponent < 0)
 1587                 abs_int_exponent = -abs_int_exponent;
 1588 
 1589             mData->mImmed.pop_back(); mData->mByteCode.pop_back();
 1590             --mStackPtr;
 1591             // ^Though the above is accounted for by the procedure
 1592             // that generates cPow, we need it for correct cFetch
 1593             // indexes in CompilePowi().
 1594 
 1595             while(sqrt_count > 0)
 1596             {
 1597                 int opcode = cSqrt;
 1598                 if(sqrt_count == 1 && int_exponent < 0)
 1599                 {
 1600                     opcode = cRSqrt;
 1601                     int_exponent = -int_exponent;
 1602                 }
 1603                 mData->mByteCode.push_back(opcode);
 1604                 --sqrt_count;
 1605             }
 1606             if((abs_int_exponent & 1) == 0)
 1607             {
 1608                 // This special rule fixes the optimization
 1609                 // shortcoming of (-x)^2 with minimal overhead.
 1610                 AddFunctionOpcode(cSqr);
 1611                 abs_int_exponent >>= 1;
 1612             }
 1613             CompilePowi(abs_int_exponent);
 1614             if(int_exponent < 0) mData->mByteCode.push_back(cInv);
 1615             ++mStackPtr; // Needed because cPow adding will assume this.
 1616             return true;
 1617         }
 1618         if(sqrt_count >= 4) break;
 1619         changed_immed += changed_immed;
 1620     }
 1621 
 1622     // When we don't know whether x >= 0, we still know that
 1623     // x^y can be safely converted into exp(y * log(x))
 1624     // when y is _not_ integer, because we know that x >= 0.
 1625     // Otherwise either expression will give a NaN.
 1626     if(/*!isInteger(original_immed) ||*/
 1627        IsNeverNegativeValueOpcode(mData->mByteCode[mData->mByteCode.size()-2]))
 1628     {
 1629         mData->mImmed.pop_back();
 1630         mData->mByteCode.pop_back();
 1631         //--mStackPtr; - accounted for by the procedure that generates cPow
 1632         AddFunctionOpcode(cLog);
 1633         AddImmedOpcode(original_immed);
 1634         //incStackPtr(); - this and the next are redundant because...
 1635         AddFunctionOpcode(cMul);
 1636         //--mStackPtr;    - ...because the cImmed was popped earlier.
 1637         AddFunctionOpcode(cExp);
 1638         return true;
 1639     }
 1640     return false;
 1641 }
 1642 
 1643 //#include "fpoptimizer/opcodename.hh"
 1644 // ^ needed only if FP_TRACE_BYTECODE_OPTIMIZATION() is used
 1645 
 1646 template<typename Value_t>
 1647 inline void FunctionParserBase<Value_t>::AddFunctionOpcode(unsigned opcode)
 1648 {
 1649 #define FP_FLOAT_VERSION 1
 1650 #define FP_COMPLEX_VERSION 0
 1651 #include "extrasrc/fp_opcode_add.inc"
 1652 #undef FP_COMPLEX_VERSION
 1653 #undef FP_FLOAT_VERSION
 1654 }
 1655 
 1656 #ifdef FP_SUPPORT_LONG_INT_TYPE
 1657 template<>
 1658 inline void FunctionParserBase<long>::AddFunctionOpcode(unsigned opcode)
 1659 {
 1660     typedef long Value_t;
 1661 #define FP_FLOAT_VERSION 0
 1662 #define FP_COMPLEX_VERSION 0
 1663 #include "extrasrc/fp_opcode_add.inc"
 1664 #undef FP_COMPLEX_VERSION
 1665 #undef FP_FLOAT_VERSION
 1666 }
 1667 #endif
 1668 
 1669 #ifdef FP_SUPPORT_GMP_INT_TYPE
 1670 template<>
 1671 inline void FunctionParserBase<GmpInt>::AddFunctionOpcode(unsigned opcode)
 1672 {
 1673     typedef GmpInt Value_t;
 1674 #define FP_FLOAT_VERSION 0
 1675 #define FP_COMPLEX_VERSION 0
 1676 #include "extrasrc/fp_opcode_add.inc"
 1677 #undef FP_COMPLEX_VERSION
 1678 #undef FP_FLOAT_VERSION
 1679 }
 1680 #endif
 1681 
 1682 #ifdef FP_SUPPORT_COMPLEX_DOUBLE_TYPE
 1683 template<>
 1684 inline void FunctionParserBase<std::complex<double> >::AddFunctionOpcode(unsigned opcode)
 1685 {
 1686     typedef std::complex<double> Value_t;
 1687 #define FP_FLOAT_VERSION 1
 1688 #define FP_COMPLEX_VERSION 1
 1689 #include "extrasrc/fp_opcode_add.inc"
 1690 #undef FP_COMPLEX_VERSION
 1691 #undef FP_FLOAT_VERSION
 1692 }
 1693 #endif
 1694 
 1695 #ifdef FP_SUPPORT_COMPLEX_FLOAT_TYPE
 1696 template<>
 1697 inline void FunctionParserBase<std::complex<float> >::AddFunctionOpcode(unsigned opcode)
 1698 {
 1699     typedef std::complex<float> Value_t;
 1700 #define FP_FLOAT_VERSION 1
 1701 #define FP_COMPLEX_VERSION 1
 1702 #include "extrasrc/fp_opcode_add.inc"
 1703 #undef FP_COMPLEX_VERSION
 1704 #undef FP_FLOAT_VERSION
 1705 }
 1706 #endif
 1707 
 1708 #ifdef FP_SUPPORT_COMPLEX_LONG_DOUBLE_TYPE
 1709 template<>
 1710 inline void FunctionParserBase<std::complex<long double> >::AddFunctionOpcode(unsigned opcode)
 1711 {
 1712     typedef std::complex<long double> Value_t;
 1713 #define FP_FLOAT_VERSION 1
 1714 #define FP_COMPLEX_VERSION 1
 1715 #include "extrasrc/fp_opcode_add.inc"
 1716 #undef FP_COMPLEX_VERSION
 1717 #undef FP_FLOAT_VERSION
 1718 }
 1719 #endif
 1720 
 1721 template<typename Value_t>
 1722 unsigned
 1723 FunctionParserBase<Value_t>::ParseIdentifier(const char* function)
 1724 {
 1725     return readIdentifier<Value_t>(function);
 1726 }
 1727 
 1728 template<typename Value_t>
 1729 std::pair<const char*, Value_t>
 1730 FunctionParserBase<Value_t>::ParseLiteral(const char* function)
 1731 {
 1732     char* endptr;
 1733 #if 0 /* Profile the hex literal parser */
 1734     if(function[0]=='0' && function[1]=='x')
 1735     {
 1736         // Parse hexadecimal literal if fp_parseLiteral didn't already
 1737         Value_t val = parseHexLiteral<Value_t>(function+2, &endptr);
 1738         if(endptr == function+2)
 1739             return std::pair<const char*,Value_t> (function, Value_t());
 1740         return std::pair<const char*, Value_t> (endptr, val);
 1741     }
 1742 #endif
 1743     Value_t val = fp_parseLiteral<Value_t>(function, &endptr);
 1744 
 1745     if(endptr == function+1 && function[0] == '0' && function[1] == 'x')
 1746     {
 1747         // Parse hexadecimal literal if fp_parseLiteral didn't already
 1748         val = parseHexLiteral<Value_t>(function+2, &endptr);
 1749         if(endptr == function+2)
 1750             return std::pair<const char*,Value_t> (function, Value_t());
 1751     }
 1752     else if(endptr == function)
 1753         return std::pair<const char*,Value_t> (function, Value_t());
 1754 
 1755     return std::pair<const char*,Value_t> (endptr, val);
 1756 }
 1757 
 1758 #ifdef FP_SUPPORT_MPFR_FLOAT_TYPE
 1759 template<>
 1760 std::pair<const char*, MpfrFloat>
 1761 FunctionParserBase<MpfrFloat>::ParseLiteral(const char* function)
 1762 {
 1763     char* endPtr;
 1764     const MpfrFloat val = MpfrFloat::parseString(function, &endPtr);
 1765     if(endPtr == function)
 1766         return std::pair<const char*,MpfrFloat> (function, MpfrFloat());
 1767     return std::pair<const char*,MpfrFloat> (endPtr, val);
 1768 }
 1769 #endif
 1770 
 1771 #ifdef FP_SUPPORT_GMP_INT_TYPE
 1772 template<>
 1773 std::pair<const char*, GmpInt>
 1774 FunctionParserBase<GmpInt>::ParseLiteral(const char* function)
 1775 {
 1776     char* endPtr;
 1777     const GmpInt val = GmpInt::parseString(function, &endPtr);
 1778     if(endPtr == function)
 1779         return std::pair<const char*,GmpInt> (function, GmpInt());
 1780     return std::pair<const char*,GmpInt> (endPtr, val);
 1781 }
 1782 #endif
 1783 
 1784 
 1785 template<typename Value_t>
 1786 inline const char*
 1787 FunctionParserBase<Value_t>::CompileLiteral(const char* function)
 1788 {
 1789     std::pair<const char*, Value_t> result = ParseLiteral(function);
 1790 
 1791     if(result.first == function)
 1792         return SetErrorType(SYNTAX_ERROR, result.first);
 1793 
 1794     AddImmedOpcode(result.second);
 1795     incStackPtr();
 1796     SkipSpace(result.first);
 1797     return result.first;
 1798 }
 1799 
 1800 template<typename Value_t>
 1801 const char* FunctionParserBase<Value_t>::CompileIf(const char* function)
 1802 {
 1803     if(*function != '(') return SetErrorType(EXPECT_PARENTH_FUNC, function);
 1804 
 1805     function = CompileExpression(function+1);
 1806     if(!function) return 0;
 1807     if(*function != ',')
 1808         return SetErrorType(noCommaError<Value_t>(*function), function);
 1809 
 1810     OPCODE opcode = cIf;
 1811     if(mData->mByteCode.back() == cNotNot) mData->mByteCode.pop_back();
 1812     if(IsNeverNegativeValueOpcode(mData->mByteCode.back()))
 1813     {
 1814         // If we know that the condition to be tested is always
 1815         // a positive value (such as when produced by "x<y"),
 1816         // we can use the faster opcode to evaluate it.
 1817         // cIf tests whether fabs(cond) >= 0.5,
 1818         // cAbsIf simply tests whether cond >= 0.5.
 1819         opcode = cAbsIf;
 1820     }
 1821 
 1822     mData->mByteCode.push_back(opcode);
 1823     const unsigned curByteCodeSize = unsigned(mData->mByteCode.size());
 1824     PushOpcodeParam<false>(0); // Jump index; to be set later
 1825     PushOpcodeParam<true> (0); // Immed jump index; to be set later
 1826 
 1827     --mStackPtr;
 1828 
 1829     function = CompileExpression(function + 1);
 1830     if(!function) return 0;
 1831     if(*function != ',')
 1832         return SetErrorType(noCommaError<Value_t>(*function), function);
 1833 
 1834     mData->mByteCode.push_back(cJump);
 1835     const unsigned curByteCodeSize2 = unsigned(mData->mByteCode.size());
 1836     const unsigned curImmedSize2 = unsigned(mData->mImmed.size());
 1837     PushOpcodeParam<false>(0); // Jump index; to be set later
 1838     PushOpcodeParam<true> (0); // Immed jump index; to be set later
 1839 
 1840     --mStackPtr;
 1841 
 1842     function = CompileExpression(function + 1);
 1843     if(!function) return 0;
 1844     if(*function != ')')
 1845         return SetErrorType(noParenthError<Value_t>(*function), function);
 1846 
 1847     PutOpcodeParamAt<true> ( mData->mByteCode.back(), unsigned(mData->mByteCode.size()-1) );
 1848     // ^Necessary for guarding against if(x,1,2)+1 being changed
 1849     //  into if(x,1,3) by fp_opcode_add.inc
 1850 
 1851     // Set jump indices
 1852     PutOpcodeParamAt<false>( curByteCodeSize2+1, curByteCodeSize );
 1853     PutOpcodeParamAt<false>( curImmedSize2,      curByteCodeSize+1 );
 1854     PutOpcodeParamAt<false>( unsigned(mData->mByteCode.size())-1, curByteCodeSize2);
 1855     PutOpcodeParamAt<false>( unsigned(mData->mImmed.size()),      curByteCodeSize2+1);
 1856 
 1857     ++function;
 1858     SkipSpace(function);
 1859     return function;
 1860 }
 1861 
 1862 template<typename Value_t>
 1863 const char* FunctionParserBase<Value_t>::CompileFunctionParams
 1864 (const char* function, unsigned requiredParams)
 1865 {
 1866     if(*function != '(') return SetErrorType(EXPECT_PARENTH_FUNC, function);
 1867 
 1868     if(requiredParams > 0)
 1869     {
 1870         const char* function_end = CompileExpression(function+1);
 1871         if(!function_end)
 1872         {
 1873             // If an error occurred, verify whether it was caused by ()
 1874             ++function;
 1875             SkipSpace(function);
 1876             if(*function == ')')
 1877                 return SetErrorType(ILL_PARAMS_AMOUNT, function);
 1878             // Not caused by (), use the error message given by CompileExpression()
 1879             return 0;
 1880         }
 1881         function = function_end;
 1882 
 1883         for(unsigned i = 1; i < requiredParams; ++i)
 1884         {
 1885             if(*function != ',')
 1886                 return SetErrorType(noCommaError<Value_t>(*function), function);
 1887 
 1888             function = CompileExpression(function+1);
 1889             if(!function) return 0;
 1890         }
 1891         // No need for incStackPtr() because each parse parameter calls it
 1892         mStackPtr -= requiredParams-1;
 1893     }
 1894     else
 1895     {
 1896         incStackPtr(); // return value of function is pushed onto the stack
 1897         ++function;
 1898         SkipSpace(function);
 1899     }
 1900 
 1901     if(*function != ')')
 1902         return SetErrorType(noParenthError<Value_t>(*function), function);
 1903     ++function;
 1904     SkipSpace(function);
 1905     return function;
 1906 }
 1907 
 1908 template<typename Value_t>
 1909 const char* FunctionParserBase<Value_t>::CompileElement(const char* function)
 1910 {
 1911     if(BeginsLiteral<Value_t>( (unsigned char) *function))
 1912         return CompileLiteral(function);
 1913 
 1914     unsigned nameLength = readIdentifier<Value_t>(function);
 1915     if(nameLength == 0)
 1916     {
 1917         // No identifier found
 1918         if(*function == '(') return CompileParenthesis(function);
 1919         if(*function == ')') return SetErrorType(MISM_PARENTH, function);
 1920         return SetErrorType(SYNTAX_ERROR, function);
 1921     }
 1922 
 1923     // Function, variable or constant
 1924     if(nameLength & 0x80000000U) // Function
 1925     {
 1926         OPCODE func_opcode = OPCODE( (nameLength >> 16) & 0x7FFF );
 1927         return CompileFunction(function + (nameLength & 0xFFFF), func_opcode);
 1928     }
 1929 
 1930     NamePtr name(function, nameLength);
 1931     const char* endPtr = function + nameLength;
 1932     SkipSpace(endPtr);
 1933 
 1934     typename NamePtrsMap<Value_t>::iterator nameIter =
 1935         mData->mNamePtrs.find(name);
 1936     if(nameIter == mData->mNamePtrs.end())
 1937     {
 1938         // Check if it's an inline variable:
 1939         for(typename Data::InlineVarNamesContainer::reverse_iterator iter =
 1940                 mData->mInlineVarNames.rbegin();
 1941             iter != mData->mInlineVarNames.rend();
 1942             ++iter)
 1943         {
 1944             if(name == iter->mName)
 1945             {
 1946                 if( iter->mFetchIndex+1 == mStackPtr)
 1947                 {
 1948                     mData->mByteCode.push_back(cDup);
 1949                 }
 1950                 else
 1951                 {
 1952                     mData->mByteCode.push_back(cFetch);
 1953                     PushOpcodeParam<true>(iter->mFetchIndex);
 1954                 }
 1955                 incStackPtr();
 1956                 return endPtr;
 1957             }
 1958         }
 1959 
 1960         return SetErrorType(UNKNOWN_IDENTIFIER, function);
 1961     }
 1962 
 1963     const NameData<Value_t>* nameData = &nameIter->second;
 1964     switch(nameData->type)
 1965     {
 1966       case NameData<Value_t>::VARIABLE: // is variable
 1967           if(unlikely(!mData->mByteCode.empty() &&
 1968                       mData->mByteCode.back() == nameData->index))
 1969               mData->mByteCode.push_back(cDup);
 1970           else
 1971               mData->mByteCode.push_back(nameData->index);
 1972           incStackPtr();
 1973           return endPtr;
 1974 
 1975       case NameData<Value_t>::CONSTANT: // is constant
 1976           AddImmedOpcode(nameData->value);
 1977           incStackPtr();
 1978           return endPtr;
 1979 
 1980       case NameData<Value_t>::UNIT: // is unit (error if appears here)
 1981           break;
 1982 
 1983       case NameData<Value_t>::FUNC_PTR: // is C++ function
 1984           function = CompileFunctionParams
 1985               (endPtr, mData->mFuncPtrs[nameData->index].mParams);
 1986           //if(!function) return 0;
 1987           mData->mByteCode.push_back(cFCall);
 1988           PushOpcodeParam<true>(nameData->index);
 1989           return function;
 1990 
 1991       case NameData<Value_t>::PARSER_PTR: // is FunctionParser
 1992           function = CompileFunctionParams
 1993               (endPtr, mData->mFuncParsers[nameData->index].mParams);
 1994           //if(!function) return 0;
 1995           mData->mByteCode.push_back(cPCall);
 1996           PushOpcodeParam<true>(nameData->index);
 1997           return function;
 1998     }
 1999 
 2000     // When it's an unit (or unrecognized type):
 2001     return SetErrorType(SYNTAX_ERROR, function);
 2002 }
 2003 
 2004 template<typename Value_t>
 2005 inline const char* FunctionParserBase<Value_t>::CompileFunction
 2006 (const char* function, unsigned func_opcode)
 2007 {
 2008     SkipSpace(function);
 2009     const FuncDefinition& funcDef = Functions[func_opcode];
 2010 
 2011     if(func_opcode == cIf) // "if" is a special case
 2012         return CompileIf(function);
 2013 
 2014     unsigned requiredParams = funcDef.params;
 2015 
 2016     function = CompileFunctionParams(function, requiredParams);
 2017     if(!function) return 0;
 2018 
 2019     if(mData->mUseDegreeConversion)
 2020     {
 2021         if(funcDef.flags & FuncDefinition::AngleIn)
 2022             AddFunctionOpcode(cRad);
 2023 
 2024         AddFunctionOpcode(func_opcode);
 2025 
 2026         if(funcDef.flags & FuncDefinition::AngleOut)
 2027             AddFunctionOpcode(cDeg);
 2028     }
 2029     else
 2030     {
 2031         AddFunctionOpcode(func_opcode);
 2032     }
 2033     return function;
 2034 }
 2035 
 2036 template<typename Value_t>
 2037 inline const char*
 2038 FunctionParserBase<Value_t>::CompileParenthesis(const char* function)
 2039 {
 2040     ++function; // Skip '('
 2041 
 2042     SkipSpace(function);
 2043     if(*function == ')') return SetErrorType(EMPTY_PARENTH, function);
 2044     function = CompileExpression(function);
 2045     if(!function) return 0;
 2046 
 2047     if(*function != ')') return SetErrorType(MISSING_PARENTH, function);
 2048     ++function; // Skip ')'
 2049 
 2050     SkipSpace(function);
 2051     return function;
 2052 }
 2053 
 2054 template<typename Value_t>
 2055 const char*
 2056 FunctionParserBase<Value_t>::CompilePossibleUnit(const char* function)
 2057 {
 2058     unsigned nameLength = readIdentifier<Value_t>(function);
 2059     if(nameLength & 0x80000000U) return function; // built-in function name
 2060     if(nameLength != 0)
 2061     {
 2062         NamePtr name(function, nameLength);
 2063 
 2064         typename NamePtrsMap<Value_t>::iterator nameIter =
 2065             mData->mNamePtrs.find(name);
 2066         if(nameIter != mData->mNamePtrs.end())
 2067         {
 2068             const NameData<Value_t>* nameData = &nameIter->second;
 2069             if(nameData->type == NameData<Value_t>::UNIT)
 2070             {
 2071                 AddImmedOpcode(nameData->value);
 2072                 incStackPtr();
 2073                 AddFunctionOpcode(cMul);
 2074                 --mStackPtr;
 2075 
 2076                 const char* endPtr = function + nameLength;
 2077                 SkipSpace(endPtr);
 2078                 return endPtr;
 2079             }
 2080         }
 2081     }
 2082 
 2083     return function;
 2084 }
 2085 
 2086 template<typename Value_t>
 2087 inline const char*
 2088 FunctionParserBase<Value_t>::CompilePow(const char* function)
 2089 {
 2090     function = CompileElement(function);
 2091     if(!function) return 0;
 2092     function = CompilePossibleUnit(function);
 2093 
 2094     if(*function == '^')
 2095     {
 2096         ++function;
 2097         SkipSpace(function);
 2098 
 2099         unsigned op = cPow;
 2100         if(mData->mByteCode.back() == cImmed)
 2101         {
 2102             if(mData->mImmed.back() == fp_const_e<Value_t>())
 2103                 { op = cExp;  mData->mByteCode.pop_back();
 2104                     mData->mImmed.pop_back(); --mStackPtr; }
 2105             else if(mData->mImmed.back() == Value_t(2))
 2106                 { op = cExp2; mData->mByteCode.pop_back();
 2107                     mData->mImmed.pop_back(); --mStackPtr; }
 2108         }
 2109 
 2110         function = CompileUnaryMinus(function);
 2111         if(!function) return 0;
 2112 
 2113         // add opcode
 2114         AddFunctionOpcode(op);
 2115 
 2116         if(op == cPow) --mStackPtr;
 2117     }
 2118     return function;
 2119 }
 2120 
 2121 /* Currently the power operator is skipped for integral types because its
 2122    usefulness with them is questionable, and in the case of GmpInt, for safety
 2123    reasons:
 2124    - With long int almost any power, except for very small ones, would
 2125      overflow the result, so the usefulness of this is rather questionable.
 2126    - With GmpInt the power operator could be easily abused to make the program
 2127      run out of memory (think of a function like "10^10^10^10^1000000").
 2128 */
 2129 #ifdef FP_SUPPORT_LONG_INT_TYPE
 2130 template<>
 2131 inline const char*
 2132 FunctionParserBase<long>::CompilePow(const char* function)
 2133 {
 2134     function = CompileElement(function);
 2135     if(!function) return 0;
 2136     return CompilePossibleUnit(function);
 2137 }
 2138 #endif
 2139 
 2140 #ifdef FP_SUPPORT_GMP_INT_TYPE
 2141 template<>
 2142 inline const char*
 2143 FunctionParserBase<GmpInt>::CompilePow(const char* function)
 2144 {
 2145     function = CompileElement(function);
 2146     if(!function) return 0;
 2147     return CompilePossibleUnit(function);
 2148 }
 2149 #endif
 2150 
 2151 template<typename Value_t>
 2152 inline const char*
 2153 FunctionParserBase<Value_t>::CompileUnaryMinus(const char* function)
 2154 {
 2155     char op = *function;
 2156     switch(op)
 2157     {
 2158         case '-':
 2159         case '!':
 2160             ++function;
 2161             SkipSpace(function);
 2162 
 2163             function = CompileUnaryMinus(function);
 2164             if(!function) return 0;
 2165 
 2166             AddFunctionOpcode(op=='-' ? cNeg : cNot);
 2167 
 2168             return function;
 2169         default: break;
 2170     }
 2171     return CompilePow(function);
 2172 }
 2173 
 2174 template<typename Value_t>
 2175 inline const char*
 2176 FunctionParserBase<Value_t>::CompileMult(const char* function)
 2177 {
 2178     function = CompileUnaryMinus(function);
 2179     if(!function) return 0;
 2180 
 2181     Value_t pending_immed(1);
 2182     #define FP_FlushImmed(do_reset) \
 2183         if(pending_immed != Value_t(1)) \
 2184         { \
 2185             unsigned op = cMul; \
 2186             if(!IsIntType<Value_t>::result && mData->mByteCode.back() == cInv) \
 2187             { \
 2188                 /* (...) cInv 5 cMul -> (...) 5 cRDiv */ \
 2189                 /*           ^               ^      | */ \
 2190                 mData->mByteCode.pop_back(); \
 2191                 op = cRDiv; \
 2192             } \
 2193             AddImmedOpcode(pending_immed); \
 2194             incStackPtr(); \
 2195             AddFunctionOpcode(op); \
 2196             --mStackPtr; \
 2197             if(do_reset) pending_immed = Value_t(1); \
 2198         }
 2199     while(true)
 2200     {
 2201         char c = *function;
 2202         if(c == '%')
 2203         {
 2204             FP_FlushImmed(true);
 2205             ++function;
 2206             SkipSpace(function);
 2207             function = CompileUnaryMinus(function);
 2208             if(!function) return 0;
 2209             AddFunctionOpcode(cMod);
 2210             --mStackPtr;
 2211             continue;
 2212         }
 2213         if(c != '*' && c != '/') break;
 2214 
 2215         bool safe_cumulation = (c == '*' || !IsIntType<Value_t>::result);
 2216         if(!safe_cumulation)
 2217         {
 2218             FP_FlushImmed(true);
 2219         }
 2220 
 2221         ++function;
 2222         SkipSpace(function);
 2223         if(mData->mByteCode.back() == cImmed
 2224         && (safe_cumulation
 2225          || mData->mImmed.back() == Value_t(1)))
 2226         {
 2227             // 5 (...) cMul --> (...)      ||| 5 cMul
 2228             // 5 (...) cDiv --> (...) cInv ||| 5 cMul
 2229             //  ^          |              ^
 2230             pending_immed *= mData->mImmed.back();
 2231             mData->mImmed.pop_back();
 2232             mData->mByteCode.pop_back();
 2233             --mStackPtr;
 2234             function = CompileUnaryMinus(function);
 2235             if(!function) return 0;
 2236             if(c == '/')
 2237                 AddFunctionOpcode(cInv);
 2238             continue;
 2239         }
 2240         if(safe_cumulation
 2241         && mData->mByteCode.back() == cMul
 2242         && mData->mByteCode[mData->mByteCode.size()-2] == cImmed)
 2243         {
 2244             // (:::) 5 cMul (...) cMul -> (:::) (...) cMul  ||| 5 cMul
 2245             // (:::) 5 cMul (...) cDiv -> (:::) (...) cDiv  ||| 5 cMul
 2246             //             ^                   ^
 2247             pending_immed *= mData->mImmed.back();
 2248             mData->mImmed.pop_back();
 2249             mData->mByteCode.pop_back();
 2250             mData->mByteCode.pop_back();
 2251         }
 2252         // cDiv is not tested here because the bytecode
 2253         // optimizer will convert this kind of cDivs into cMuls.
 2254         bool lhs_inverted = false;
 2255         if(!IsIntType<Value_t>::result && c == '*'
 2256         && mData->mByteCode.back() == cInv)
 2257         {
 2258             // (:::) cInv (...) cMul -> (:::) (...) cRDiv
 2259             // (:::) cInv (...) cDiv -> (:::) (...) cMul cInv
 2260             //           ^                   ^            |
 2261             mData->mByteCode.pop_back();
 2262             lhs_inverted = true;
 2263         }
 2264         function = CompileUnaryMinus(function);
 2265         if(!function) return 0;
 2266         if(safe_cumulation
 2267         && mData->mByteCode.back() == cMul
 2268         && mData->mByteCode[mData->mByteCode.size()-2] == cImmed)
 2269         {
 2270             // (:::) (...) 5 cMul cMul -> (:::) (...) cMul  |||  5 Mul
 2271             // (:::) (...) 5 cMul cDiv -> (:::) (...) cDiv  ||| /5 Mul
 2272             //                   ^                        ^
 2273             if(c == '*')
 2274                 pending_immed *= mData->mImmed.back();
 2275             else
 2276                 pending_immed /= mData->mImmed.back();
 2277             mData->mImmed.pop_back();
 2278             mData->mByteCode.pop_back();
 2279             mData->mByteCode.pop_back();
 2280         }
 2281         else
 2282         if(safe_cumulation
 2283         && mData->mByteCode.back() == cRDiv
 2284         && mData->mByteCode[mData->mByteCode.size()-2] == cImmed)
 2285         {
 2286             // (:::) (...) 5 cRDiv cMul -> (:::) (...) cDiv  |||  5 cMul
 2287             // (:::) (...) 5 cRDiv cDiv -> (:::) (...) cMul  ||| /5 cMul
 2288             //                    ^                   ^
 2289             if(c == '*')
 2290                 { c = '/'; pending_immed *= mData->mImmed.back(); }
 2291             else
 2292                 { c = '*'; pending_immed /= mData->mImmed.back(); }
 2293             mData->mImmed.pop_back();
 2294             mData->mByteCode.pop_back();
 2295             mData->mByteCode.pop_back();
 2296         }
 2297         if(!lhs_inverted) // if (/x/y) was changed to /(x*y), add missing cInv
 2298         {
 2299             AddFunctionOpcode(c == '*' ? cMul : cDiv);
 2300             --mStackPtr;
 2301         }
 2302         else if(c == '*') // (/x)*y -> rdiv(x,y)
 2303         {
 2304             AddFunctionOpcode(cRDiv);
 2305             --mStackPtr;
 2306         }
 2307         else // (/x)/y -> /(x*y)
 2308         {
 2309             AddFunctionOpcode(cMul);
 2310             --mStackPtr;
 2311             AddFunctionOpcode(cInv);
 2312         }
 2313     }
 2314     FP_FlushImmed(false);
 2315     #undef FP_FlushImmed
 2316     return function;
 2317 }
 2318 
 2319 template<typename Value_t>
 2320 inline const char*
 2321 FunctionParserBase<Value_t>::CompileAddition(const char* function)
 2322 {
 2323     function = CompileMult(function);
 2324     if(!function) return 0;
 2325 
 2326     Value_t pending_immed(0);
 2327     #define FP_FlushImmed(do_reset) \
 2328         if(pending_immed != Value_t(0)) \
 2329         { \
 2330             unsigned op = cAdd; \
 2331             if(mData->mByteCode.back() == cNeg) \
 2332             { \
 2333                 /* (...) cNeg 5 cAdd -> (...) 5 cRSub */ \
 2334                 /*           ^               ^      | */ \
 2335                 mData->mByteCode.pop_back(); \
 2336                 op = cRSub; \
 2337             } \
 2338             AddImmedOpcode(pending_immed); \
 2339             incStackPtr(); \
 2340             AddFunctionOpcode(op); \
 2341             --mStackPtr; \
 2342             if(do_reset) pending_immed = Value_t(0); \
 2343         }
 2344     while(true)
 2345     {
 2346         char c = *function;
 2347         if(c != '+' && c != '-') break;
 2348         ++function;
 2349         SkipSpace(function);
 2350         if(mData->mByteCode.back() == cImmed)
 2351         {
 2352             // 5 (...) cAdd --> (...)      ||| 5 cAdd
 2353             // 5 (...) cSub --> (...) cNeg ||| 5 cAdd
 2354             //  ^          |              ^
 2355             pending_immed += mData->mImmed.back();
 2356             mData->mImmed.pop_back();
 2357             mData->mByteCode.pop_back();
 2358             --mStackPtr;
 2359             function = CompileMult(function);
 2360             if(!function) return 0;
 2361             if(c == '-')
 2362                 AddFunctionOpcode(cNeg);
 2363             continue;
 2364         }
 2365         if(mData->mByteCode.back() == cAdd
 2366         && mData->mByteCode[mData->mByteCode.size()-2] == cImmed)
 2367         {
 2368             // (:::) 5 cAdd (...) cAdd -> (:::) (...) cAdd  ||| 5 cAdd
 2369             // (:::) 5 cAdd (...) cSub -> (:::) (...) cSub  ||| 5 cAdd
 2370             //             ^                   ^
 2371             pending_immed += mData->mImmed.back();
 2372             mData->mImmed.pop_back();
 2373             mData->mByteCode.pop_back();
 2374             mData->mByteCode.pop_back();
 2375         }
 2376         // cSub is not tested here because the bytecode
 2377         // optimizer will convert this kind of cSubs into cAdds.
 2378         bool lhs_negated = false;
 2379         if(mData->mByteCode.back() == cNeg)
 2380         {
 2381             // (:::) cNeg (...) cAdd -> (:::) (...) cRSub
 2382             // (:::) cNeg (...) cSub -> (:::) (...) cAdd cNeg
 2383             //           ^                   ^            |
 2384             mData->mByteCode.pop_back();
 2385             lhs_negated = true;
 2386         }
 2387         function = CompileMult(function);
 2388         if(!function) return 0;
 2389         if(mData->mByteCode.back() == cAdd
 2390         && mData->mByteCode[mData->mByteCode.size()-2] == cImmed)
 2391         {
 2392             // (:::) (...) 5 cAdd cAdd -> (:::) (...) cAdd  |||  5 Add
 2393             // (:::) (...) 5 cAdd cSub -> (:::) (...) cSub  ||| -5 Add
 2394             //                   ^                        ^
 2395             if(c == '+')
 2396                 pending_immed += mData->mImmed.back();
 2397             else
 2398                 pending_immed -= mData->mImmed.back();
 2399             mData->mImmed.pop_back();
 2400             mData->mByteCode.pop_back();
 2401             mData->mByteCode.pop_back();
 2402         }
 2403         else
 2404         if(mData->mByteCode.back() == cRSub
 2405         && mData->mByteCode[mData->mByteCode.size()-2] == cImmed)
 2406         {
 2407             // (:::) (...) 5 cRSub cAdd -> (:::) (...) cSub  |||  5 cAdd
 2408             // (:::) (...) 5 cRSub cSub -> (:::) (...) cAdd  ||| -5 cAdd
 2409             //                    ^                   ^
 2410             if(c == '+')
 2411                 { c = '-'; pending_immed += mData->mImmed.back(); }
 2412             else
 2413                 { c = '+'; pending_immed -= mData->mImmed.back(); }
 2414             mData->mImmed.pop_back();
 2415             mData->mByteCode.pop_back();
 2416             mData->mByteCode.pop_back();
 2417         }
 2418         if(!lhs_negated) // if (-x-y) was changed to -(x+y), add missing cNeg
 2419         {
 2420             AddFunctionOpcode(c == '+' ? cAdd : cSub);
 2421             --mStackPtr;
 2422         }
 2423         else if(c == '+') // (-x)+y -> rsub(x,y)
 2424         {
 2425             AddFunctionOpcode(cRSub);
 2426             --mStackPtr;
 2427         }
 2428         else // (-x)-y -> -(x+y)
 2429         {
 2430             AddFunctionOpcode(cAdd);
 2431             --mStackPtr;
 2432             AddFunctionOpcode(cNeg);
 2433         }
 2434     }
 2435     FP_FlushImmed(false);
 2436     #undef FP_FlushImmed
 2437     return function;
 2438 }
 2439 
 2440 template<typename Value_t>
 2441 inline const char*
 2442 FunctionParserBase<Value_t>::CompileComparison(const char* function)
 2443 {
 2444     unsigned op=0;
 2445     while(true)
 2446     {
 2447         function = CompileAddition(function);
 2448         if(!function) return 0;
 2449 
 2450         if(op)
 2451         {
 2452             AddFunctionOpcode(op);
 2453             --mStackPtr;
 2454         }
 2455         switch(*function)
 2456         {
 2457           case '=':
 2458               ++function; op = cEqual; break;
 2459           case '!':
 2460               if(function[1] == '=')
 2461               { function += 2; op = cNEqual; break; }
 2462               // If '=' does not follow '!', a syntax error will
 2463               // be generated at the outermost parsing level
 2464               return function;
 2465           case '<':
 2466               if(function[1] == '=')
 2467               { function += 2; op = cLessOrEq; break; }
 2468               ++function; op = cLess; break;
 2469           case '>':
 2470               if(function[1] == '=')
 2471               { function += 2; op = cGreaterOrEq; break; }
 2472               ++function; op = cGreater; break;
 2473           default: return function;
 2474         }
 2475         SkipSpace(function);
 2476     }
 2477     return function;
 2478 }
 2479 
 2480 template<typename Value_t>
 2481 inline const char* FunctionParserBase<Value_t>::CompileAnd(const char* function)
 2482 {
 2483     std::size_t param0end=0;
 2484     while(true)
 2485     {
 2486         function = CompileComparison(function);
 2487         if(!function) return 0;
 2488 
 2489         if(param0end)
 2490         {
 2491             if(mData->mByteCode.back() == cNotNot) mData->mByteCode.pop_back();
 2492 
 2493             AddFunctionOpcode(cAnd);
 2494             --mStackPtr;
 2495         }
 2496         if(*function != '&') break;
 2497         ++function;
 2498         SkipSpace(function);
 2499         param0end = mData->mByteCode.size();
 2500     }
 2501     return function;
 2502 }
 2503 
 2504 template<typename Value_t>
 2505 const char* FunctionParserBase<Value_t>::CompileExpression(const char* function)
 2506 {
 2507     std::size_t param0end=0;
 2508     while(true)
 2509     {
 2510         SkipSpace(function);
 2511         function = CompileAnd(function);
 2512         if(!function) return 0;
 2513 
 2514         if(param0end)
 2515         {
 2516             if(mData->mByteCode.back() == cNotNot) mData->mByteCode.pop_back();
 2517 
 2518             AddFunctionOpcode(cOr);
 2519             --mStackPtr;
 2520         }
 2521         if(*function != '|') break;
 2522         ++function;
 2523         param0end = mData->mByteCode.size();
 2524     }
 2525     return function;
 2526 }
 2527 
 2528 template<typename Value_t>
 2529 const char* FunctionParserBase<Value_t>::Compile(const char* function)
 2530 {
 2531     while(true)
 2532     {
 2533         // Check if an identifier appears as first token:
 2534         SkipSpace(function);
 2535         unsigned nameLength = readIdentifier<Value_t>(function);
 2536         if(nameLength > 0 && !(nameLength & 0x80000000U))
 2537         {
 2538             typename Data::InlineVariable inlineVar =
 2539                 { NamePtr(function, nameLength), 0 };
 2540 
 2541             // Check if it's an unknown identifier:
 2542             typename NamePtrsMap<Value_t>::iterator nameIter =
 2543                 mData->mNamePtrs.find(inlineVar.mName);
 2544             if(nameIter == mData->mNamePtrs.end())
 2545             {
 2546                 const char* function2 = function + nameLength;
 2547                 SkipSpace(function2);
 2548 
 2549                 // Check if ":=" follows the unknown identifier:
 2550                 if(function2[0] == ':' && function2[1] == '=')
 2551                 {
 2552                     // Parse the expression that follows and create the
 2553                     // inline variable:
 2554                     function2 = CompileExpression(function2 + 2);
 2555                     if(!function2) return 0;
 2556                     if(*function2 != ';') return function2;
 2557 
 2558                     inlineVar.mFetchIndex = mStackPtr - 1;
 2559                     mData->mInlineVarNames.push_back(inlineVar);
 2560 
 2561                     // Continue with the expression after the ';':
 2562                     function = function2 + 1;
 2563                     continue;
 2564                 }
 2565             }
 2566         }
 2567         break;
 2568     }
 2569 
 2570     return CompileExpression(function);
 2571 }
 2572 
 2573 template<typename Value_t> template<bool PutFlag>
 2574 inline void FunctionParserBase<Value_t>::PushOpcodeParam
 2575     (unsigned value)
 2576 {
 2577     mData->mByteCode.push_back(value | (PutFlag ? FP_ParamGuardMask : 0u));
 2578     if(PutFlag) mData->mHasByteCodeFlags = true;
 2579 }
 2580 
 2581 template<typename Value_t> template<bool PutFlag>
 2582 inline void FunctionParserBase<Value_t>::PutOpcodeParamAt
 2583     (unsigned value, unsigned offset)
 2584 {
 2585     mData->mByteCode[offset] = value | (PutFlag ? FP_ParamGuardMask : 0u);
 2586     if(PutFlag) mData->mHasByteCodeFlags = true;
 2587 }
 2588 
 2589 //===========================================================================
 2590 // Function evaluation
 2591 //===========================================================================
 2592 template<typename Value_t>
 2593 Value_t FunctionParserBase<Value_t>::Eval(const Value_t* Vars)
 2594 {
 2595     if(mData->mParseErrorType != FP_NO_ERROR) return Value_t(0);
 2596 
 2597     const unsigned* const byteCode = &(mData->mByteCode[0]);
 2598     const Value_t* const immed = mData->mImmed.empty() ? 0 : &(mData->mImmed[0]);
 2599     const unsigned byteCodeSize = unsigned(mData->mByteCode.size());
 2600     unsigned IP, DP=0;
 2601     int SP=-1;
 2602     int stt =0;
 2603     std::vector<Value_t>& StackSave = mData->mStackSave;
 2604     int VarSize = StackSave.size();
 2605 
 2606 #ifdef FP_USE_THREAD_SAFE_EVAL
 2607     /* If Eval() may be called by multiple threads simultaneously,
 2608      * then Eval() must allocate its own stack.
 2609      */
 2610 #ifdef FP_USE_THREAD_SAFE_EVAL_WITH_ALLOCA
 2611     /* alloca() allocates room from the hardware stack.
 2612      * It is automatically freed when the function returns.
 2613      */
 2614     Value_t* const Stack = (Value_t*)alloca(mData->mStackSize*sizeof(Value_t));
 2615 #else
 2616     /* Allocate from the heap. Ensure that it is freed
 2617      * automatically no matter which exit path is taken.
 2618      */
 2619     struct AutoDealloc
 2620     {
 2621         Value_t* ptr;
 2622         ~AutoDealloc() { delete[] ptr; }
 2623     } AutoDeallocStack = { new Value_t[mData->mStackSize] };
 2624     Value_t*& Stack = AutoDeallocStack.ptr;
 2625 #endif
 2626 #else
 2627     /* No thread safety, so use a global stack. */
 2628     std::vector<Value_t>& Stack = mData->mStack;
 2629 #endif
 2630 
 2631     for(IP=0; IP<byteCodeSize; ++IP)
 2632     {
 2633         switch(byteCode[IP])
 2634         {
 2635 // Functions:
 2636           case   cAbs: Stack[SP] = fp_abs(Stack[SP]); break;
 2637 
 2638           case  cAcos:
 2639                 if(IsComplexType<Value_t>::result == false
 2640                 && (Stack[SP] < Value_t(-1) || Stack[SP] > Value_t(1)))
 2641                 {
 2642                   mData->mEvalErrorType=4;
 2643                   return Value_t(0);
 2644                 }
 2645                 Stack[SP] = fp_acos(Stack[SP]);
 2646                 break;
 2647 
 2648           case cAcosh:
 2649             if(IsComplexType<Value_t>::result == false
 2650             && Stack[SP] < Value_t(1))
 2651             {
 2652               mData->mEvalErrorType=4;
 2653               return Value_t(0);
 2654             }
 2655             Stack[SP] = fp_acosh(Stack[SP]);
 2656             break;
 2657 
 2658           case  cAsin:
 2659             if(IsComplexType<Value_t>::result == false
 2660             && (Stack[SP] < Value_t(-1) || Stack[SP] > Value_t(1)))
 2661             {
 2662               mData->mEvalErrorType=4;
 2663               return Value_t(0);
 2664             }
 2665             Stack[SP] = fp_asin(Stack[SP]);
 2666             break;
 2667 
 2668           case cAsinh:
 2669             Stack[SP] = fp_asinh(Stack[SP]);
 2670             break;
 2671 
 2672           case  cAtan: Stack[SP] = fp_atan(Stack[SP]); break;
 2673 
 2674           case cAtan2: Stack[SP-1] = fp_atan2(Stack[SP-1], Stack[SP]);
 2675                        --SP; break;
 2676 
 2677           case cAtanh:
 2678               if(IsComplexType<Value_t>::result
 2679               ?  (Stack[SP] == Value_t(-1) || Stack[SP] == Value_t(1))
 2680               :  (Stack[SP] <= Value_t(-1) || Stack[SP] >= Value_t(1)))
 2681               { mData->mEvalErrorType=4; return Value_t(0); }
 2682               Stack[SP] = fp_atanh(Stack[SP]); break;
 2683 
 2684           case  cCbrt: Stack[SP] = fp_cbrt(Stack[SP]); break;
 2685 
 2686           case  cCeil: Stack[SP] = fp_ceil(Stack[SP]); break;
 2687 
 2688           case   cCos: Stack[SP] = fp_cos(Stack[SP]); break;
 2689 
 2690           case  cCosh: Stack[SP] = fp_cosh(Stack[SP]); break;
 2691 
 2692           case   cCot:
 2693               {
 2694                   const Value_t t = fp_tan(Stack[SP]);
 2695                   if(t == Value_t(0))
 2696                   { mData->mEvalErrorType=1; return Value_t(0); }
 2697                   Stack[SP] = Value_t(1)/t; break;
 2698               }
 2699 
 2700           case   cCsc:
 2701               {
 2702                   const Value_t s = fp_sin(Stack[SP]);
 2703                   if(s == Value_t(0))
 2704                   { mData->mEvalErrorType=1; return Value_t(0); }
 2705                   Stack[SP] = Value_t(1)/s; break;
 2706               }
 2707 
 2708 
 2709           case   cExp: Stack[SP] = fp_exp(Stack[SP]); break;
 2710 
 2711           case   cExp2: Stack[SP] = fp_exp2(Stack[SP]); break;
 2712 
 2713           case cFloor: Stack[SP] = fp_floor(Stack[SP]); break;
 2714 
 2715           case cHypot:
 2716               Stack[SP-1] = fp_hypot(Stack[SP-1], Stack[SP]);
 2717               --SP; break;
 2718 
 2719           case    cIf:
 2720                   if(fp_truth(Stack[SP--]))
 2721                       IP += 2;
 2722                   else
 2723                   {
 2724                       const unsigned* buf = &byteCode[IP+1];
 2725                       IP = buf[0];
 2726                       DP = buf[1];
 2727                   }
 2728                   break;
 2729 
 2730           case   cInt: Stack[SP] = fp_int(Stack[SP]); break;
 2731 
 2732           case   cLog:
 2733               if(IsComplexType<Value_t>::result
 2734                ?   Stack[SP] == Value_t(0)
 2735                :   !(Stack[SP] > Value_t(0)))
 2736               { mData->mEvalErrorType=3; return Value_t(0); }
 2737               Stack[SP] = fp_log(Stack[SP]); break;
 2738 
 2739           case cLog10:
 2740               if(IsComplexType<Value_t>::result
 2741                ?   Stack[SP] == Value_t(0)
 2742                :   !(Stack[SP] > Value_t(0)))
 2743               { mData->mEvalErrorType=3; return Value_t(0); }
 2744               Stack[SP] = fp_log10(Stack[SP]);
 2745               break;
 2746 
 2747           case  cLog2:
 2748               if(IsComplexType<Value_t>::result
 2749                ?   Stack[SP] == Value_t(0)
 2750                :   !(Stack[SP] > Value_t(0)))
 2751               { mData->mEvalErrorType=3; return Value_t(0); }
 2752               Stack[SP] = fp_log2(Stack[SP]);
 2753               break;
 2754 
 2755           case   cMax: Stack[SP-1] = fp_max(Stack[SP-1], Stack[SP]);
 2756                        --SP; break;
 2757 
 2758           case   cMin: Stack[SP-1] = fp_min(Stack[SP-1], Stack[SP]);
 2759                        --SP; break;
 2760 
 2761           case   cPow:
 2762               // x:Negative ^ y:NonInteger is failure,
 2763               // except when the reciprocal of y forms an integer
 2764               /*if(IsComplexType<Value_t>::result == false
 2765               && Stack[SP-1] < Value_t(0) &&
 2766                  !isInteger(Stack[SP]) &&
 2767                  !isInteger(1.0 / Stack[SP]))
 2768               { mEvalErrorType=3; return Value_t(0); }*/
 2769               // x:0 ^ y:negative is failure
 2770               if(Stack[SP-1] == Value_t(0) && Stack[SP] < Value_t(0))
 2771               {
 2772                   mData->mEvalErrorType=3;
 2773                   return Value_t(0);
 2774               }
 2775               Stack[SP-1] = fp_pow(Stack[SP-1], Stack[SP]);
 2776               --SP;
 2777             break;
 2778 
 2779           case cPsh:
 2780             if((stt= Stack[SP-1]) >=VarSize || stt <0)
 2781             {
 2782                 mData->mEvalErrorType=VAR_OVERFLOW;
 2783                 return Value_t(VAR_OVERFLOW);
 2784             }
 2785             StackSave[stt] = Stack[SP];
 2786             Stack[SP-1]  = 1.0;
 2787             --SP;
 2788             break;
 2789 
 2790           case cCsd:
 2791             if((stt= Stack[SP]) >=VarSize || stt <0)
 2792             {
 2793                 mData->mEvalErrorType=VAR_OVERFLOW;
 2794                 return Value_t(VAR_OVERFLOW);
 2795             }
 2796             Stack[SP] = StackSave[stt];
 2797             break;
 2798 
 2799           case  cTrunc:
 2800             Stack[SP] = fp_trunc(Stack[SP]);
 2801             break;
 2802 
 2803           case   cSec:
 2804           {
 2805             const Value_t c = fp_cos(Stack[SP]);
 2806             if(c == Value_t(0))
 2807             {
 2808                 mData->mEvalErrorType=1;
 2809                 return Value_t(0);
 2810             }
 2811             Stack[SP] = Value_t(1)/c;
 2812             break;
 2813           }
 2814 
 2815           case   cSin:
 2816             Stack[SP] = fp_sin(Stack[SP]);
 2817             break;
 2818 
 2819           case  cSinh:
 2820             Stack[SP] = fp_sinh(Stack[SP]);
 2821             break;
 2822 
 2823           case  cSqrt:
 2824             if(IsComplexType<Value_t>::result == false &&
 2825              Stack[SP] < Value_t(0))
 2826             {
 2827               mData->mEvalErrorType=2;
 2828               return Value_t(0);
 2829             }
 2830             Stack[SP] = fp_sqrt(Stack[SP]);
 2831             break;
 2832 
 2833           case   cTan:
 2834             Stack[SP] = fp_tan(Stack[SP]);
 2835             break;
 2836 
 2837           case  cTanh:
 2838             Stack[SP] = fp_tanh(Stack[SP]);
 2839             break;
 2840 
 2841 // Misc:
 2842           case cImmed:
 2843             Stack[++SP] = immed[DP++];
 2844             break;
 2845 
 2846           case  cJump:
 2847               {
 2848                   const unsigned* buf = &byteCode[IP+1];
 2849                   IP = buf[0];
 2850                   DP = buf[1];
 2851                   break;
 2852               }
 2853 
 2854 // Operators:
 2855           case   cNeg:
 2856             Stack[SP] = -Stack[SP];
 2857             break;
 2858           case   cAdd:
 2859             Stack[SP-1] += Stack[SP];
 2860             --SP;
 2861             break;
 2862           case   cSub: Stack[SP-1] -= Stack[SP]; --SP; break;
 2863           case   cMul: Stack[SP-1] *= Stack[SP]; --SP; break;
 2864           case   cDiv:
 2865               if(Stack[SP] == Value_t(0))
 2866               { mData->mEvalErrorType=1; return Value_t(0); }
 2867               Stack[SP-1] /= Stack[SP]; --SP; break;
 2868           case   cMod:
 2869               if(Stack[SP] == Value_t(0))
 2870               { mData->mEvalErrorType=1; return Value_t(0); }
 2871               Stack[SP-1] = fp_mod(Stack[SP-1], Stack[SP]);
 2872               --SP; break;
 2873           case cEqual:
 2874               Stack[SP-1] = fp_equal(Stack[SP-1], Stack[SP]);
 2875               --SP; break;
 2876           case cNEqual:
 2877               Stack[SP-1] = fp_nequal(Stack[SP-1], Stack[SP]);
 2878               --SP; break;
 2879           case  cLess:
 2880               Stack[SP-1] = fp_less(Stack[SP-1], Stack[SP]);
 2881               --SP; break;
 2882           case  cLessOrEq:
 2883               Stack[SP-1] = fp_lessOrEq(Stack[SP-1], Stack[SP]);
 2884               --SP; break;
 2885 
 2886           case cGreater:
 2887               Stack[SP-1] = fp_less(Stack[SP], Stack[SP-1]);
 2888               --SP; break;
 2889 
 2890           case cGreaterOrEq:
 2891               Stack[SP-1] = fp_lessOrEq(Stack[SP], Stack[SP-1]);
 2892               --SP; break;
 2893           case   cNot: Stack[SP] = fp_not(Stack[SP]); break;
 2894           case cNotNot: Stack[SP] = fp_notNot(Stack[SP]); break;
 2895           case   cAnd:
 2896               Stack[SP-1] = fp_and(Stack[SP-1], Stack[SP]);
 2897               --SP; break;
 2898           case    cOr:
 2899               Stack[SP-1] = fp_or(Stack[SP-1], Stack[SP]);
 2900               --SP; break;
 2901 // Degrees-radians conversion:
 2902           case   cDeg: Stack[SP] = RadiansToDegrees(Stack[SP]); break;
 2903           case   cRad: Stack[SP] = DegreesToRadians(Stack[SP]); break;
 2904 
 2905 // User-defined function calls:
 2906           case cFCall:
 2907               {
 2908                   const unsigned index = byteCode[++IP];
 2909                   const unsigned params = mData->mFuncPtrs[index].mParams;
 2910                   const Value_t retVal =
 2911                       mData->mFuncPtrs[index].mRawFuncPtr ?
 2912                       mData->mFuncPtrs[index].mRawFuncPtr(&Stack[SP-params+1]) :
 2913                       mData->mFuncPtrs[index].mFuncWrapperPtr->callFunction
 2914                       (&Stack[SP-params+1]);
 2915                   SP -= int(params)-1;
 2916                   Stack[SP] = retVal;
 2917                   break;
 2918               }
 2919           case cPCall:
 2920               {
 2921                   unsigned index = byteCode[++IP];
 2922                   unsigned params = mData->mFuncParsers[index].mParams;
 2923                   Value_t retVal =
 2924                       mData->mFuncParsers[index].mParserPtr->Eval
 2925                       (&Stack[SP-params+1]);
 2926                   SP -= int(params)-1;
 2927                   Stack[SP] = retVal;
 2928                   const int error =
 2929                       mData->mFuncParsers[index].mParserPtr->EvalError();
 2930                   if(error)
 2931                   {
 2932                       mData->mEvalErrorType = error;
 2933                       return 0;
 2934                   }
 2935                   break;
 2936               }
 2937           case   cFetch:
 2938               {
 2939                   unsigned stackOffs = byteCode[++IP];
 2940                   Stack[SP+1] = Stack[stackOffs]; ++SP;
 2941                   break;
 2942               }
 2943 #ifdef FP_SUPPORT_OPTIMIZER
 2944           case   cPopNMov:
 2945               {
 2946                   unsigned stackOffs_target = byteCode[++IP];
 2947                   unsigned stackOffs_source = byteCode[++IP];
 2948                   Stack[stackOffs_target] = Stack[stackOffs_source];
 2949                   SP = stackOffs_target;
 2950                   break;
 2951               }
 2952 
 2953           case  cLog2by:
 2954               if(IsComplexType<Value_t>::result
 2955                ?   Stack[SP-1] == Value_t(0)
 2956                :   !(Stack[SP-1] > Value_t(0)))
 2957               { mData->mEvalErrorType=3; return Value_t(0); }
 2958               Stack[SP-1] = fp_log2(Stack[SP-1]) * Stack[SP];
 2959               --SP;
 2960               break;
 2961 
 2962           case cNop: break;
 2963 #endif // FP_SUPPORT_OPTIMIZER
 2964 
 2965           case cSinCos:
 2966               fp_sinCos(Stack[SP], Stack[SP+1], Stack[SP]);
 2967               ++SP;
 2968               break;
 2969           case cSinhCosh:
 2970               fp_sinhCosh(Stack[SP], Stack[SP+1], Stack[SP]);
 2971               ++SP;
 2972               break;
 2973 
 2974           case cAbsNot:
 2975               Stack[SP] = fp_absNot(Stack[SP]); break;
 2976           case cAbsNotNot:
 2977               Stack[SP] = fp_absNotNot(Stack[SP]); break;
 2978           case cAbsAnd:
 2979               Stack[SP-1] = fp_absAnd(Stack[SP-1], Stack[SP]);
 2980               --SP; break;
 2981           case cAbsOr:
 2982               Stack[SP-1] = fp_absOr(Stack[SP-1], Stack[SP]);
 2983               --SP; break;
 2984           case cAbsIf:
 2985               if(fp_absTruth(Stack[SP--]))
 2986                   IP += 2;
 2987               else
 2988               {
 2989                   const unsigned* buf = &byteCode[IP+1];
 2990                   IP = buf[0];
 2991                   DP = buf[1];
 2992               }
 2993               break;
 2994           case   cDup: Stack[SP+1] = Stack[SP]; ++SP; break;
 2995           case   cInv:
 2996               if(Stack[SP] == Value_t(0))
 2997               { mData->mEvalErrorType=1; return Value_t(0); }
 2998               Stack[SP] = Value_t(1)/Stack[SP];
 2999               break;
 3000           case   cSqr:
 3001               Stack[SP] = Stack[SP]*Stack[SP];
 3002               break;
 3003           case   cRDiv:
 3004               if(Stack[SP-1] == Value_t(0))
 3005               { mData->mEvalErrorType=1; return Value_t(0); }
 3006               Stack[SP-1] = Stack[SP] / Stack[SP-1]; --SP; break;
 3007           case   cRSub: Stack[SP-1] = Stack[SP] - Stack[SP-1]; --SP; break;
 3008           case   cRSqrt:
 3009               if(Stack[SP] == Value_t(0))
 3010               { mData->mEvalErrorType=1; return Value_t(0); }
 3011               Stack[SP] = Value_t(1) / fp_sqrt(Stack[SP]); break;
 3012 #ifdef FP_SUPPORT_COMPLEX_NUMBERS
 3013           case   cReal: Stack[SP] = fp_real(Stack[SP]); break;
 3014           case   cImag: Stack[SP] = fp_imag(Stack[SP]); break;
 3015           case   cArg:  Stack[SP] = fp_arg(Stack[SP]); break;
 3016           case   cConj: Stack[SP] = fp_conj(Stack[SP]); break;
 3017           case   cPolar:
 3018               Stack[SP-1] = fp_polar(Stack[SP-1], Stack[SP]);
 3019               --SP; break;
 3020 #endif
 3021 
 3022 // Variables:
 3023           default:
 3024               Stack[++SP] = Vars[byteCode[IP]-VarBegin];
 3025         }
 3026     }
 3027     mData->mEvalErrorType=0;
 3028     return Stack[SP];
 3029 }
 3030 
 3031 template<typename Value_t>
 3032 void FunctionParserBase<Value_t>::AllocateStackMemory(unsigned int nbStack, int nbvar)
 3033 {
 3034     mData->mStacki.resize(nbStack*(mData->mStack).size());
 3035     mData->mStackSave.resize(nbvar*nbStack); // Should be done earlier
 3036 }
 3037 
 3038 template<typename Value_t>
 3039 Value_t FunctionParserBase<Value_t>::Eval2(const Value_t* Vars, unsigned NbVar, double* results, unsigned NbStack, int SP)
 3040 {
 3041     if(mData->mParseErrorType != FP_NO_ERROR) return Value_t(0);
 3042 
 3043     const unsigned* const byteCode = &(mData->mByteCode[0]);
 3044     const Value_t* const immed = mData->mImmed.empty() ? 0 : &(mData->mImmed[0]);
 3045     const unsigned byteCodeSize = unsigned(mData->mByteCode.size());
 3046     unsigned IP, DP=0, Nbval;
 3047     /* No thread safety, so use a global stack. */
 3048     std::vector<Value_t>& Stack  = mData->mStack;
 3049     std::vector<Value_t>& Stacki = mData->mStacki;
 3050     std::vector<Value_t>& StackSave = mData->mStackSave;
 3051     int Size = Stack.size();
 3052     int VarSize = StackSave.size();
 3053     int stt;
 3054 
 3055     for(IP=0; IP<byteCodeSize; ++IP)
 3056     {
 3057         switch(byteCode[IP])
 3058         {
 3059 // Functions:
 3060           case   cCos:
 3061             for(Nbval=0; Nbval<NbStack; Nbval++)
 3062               Stacki[Nbval*Size+SP] = fp_cos(Stacki[Nbval*Size+SP]);
 3063             break;
 3064           case   cAbs:
 3065             for(Nbval=0; Nbval<NbStack; Nbval++)
 3066               Stacki[Nbval*Size+SP] = fp_abs(Stacki[Nbval*Size+SP]);
 3067             break;
 3068           case  cAcos:
 3069             for(Nbval=0; Nbval<NbStack; Nbval++)
 3070             {
 3071               if(IsComplexType<Value_t>::result == false
 3072               && (Stacki[Nbval*Size+SP] < Value_t(-1) || Stacki[Nbval*Size+SP] > Value_t(1)))
 3073               { mData->mEvalErrorType=4; return Value_t(0); }
 3074                   Stacki[Nbval*Size+SP] = fp_acos(Stacki[Nbval*Size+SP]);
 3075             }
 3076             break;
 3077 
 3078           case cAcosh:
 3079             for(Nbval=0; Nbval<NbStack; Nbval++)
 3080             {
 3081               if(IsComplexType<Value_t>::result == false
 3082               && Stacki[Nbval*Size+SP] < Value_t(1))
 3083               { mData->mEvalErrorType=4; return Value_t(0); }
 3084                   Stacki[Nbval*Size+SP] = fp_acosh(Stacki[Nbval*Size+SP]);
 3085             }
 3086             break;
 3087 
 3088           case  cAsin:
 3089             for(Nbval=0; Nbval<NbStack; Nbval++)
 3090             {
 3091               if(IsComplexType<Value_t>::result == false
 3092               && (Stacki[Nbval*Size+SP] < Value_t(-1) || Stacki[Nbval*Size+SP] > Value_t(1)))
 3093               { mData->mEvalErrorType=4; return Value_t(0); }
 3094                   Stacki[Nbval*Size+SP] = fp_asin(Stacki[Nbval*Size+SP]);
 3095             }
 3096             break;
 3097 
 3098           case cAsinh:
 3099             for(Nbval=0; Nbval<NbStack; Nbval++)
 3100                 Stacki[Nbval*Size+SP] = fp_asinh(Stacki[Nbval*Size+SP]);
 3101             break;
 3102 
 3103           case  cAtan:
 3104             for(Nbval=0; Nbval<NbStack; Nbval++)
 3105                 Stacki[Nbval*Size+SP] = fp_atan(Stacki[Nbval*Size+SP]);
 3106             break;
 3107 
 3108           case cAtan2:
 3109             for(Nbval=0; Nbval<NbStack; Nbval++)
 3110                 Stacki[Nbval*Size+SP-1] = fp_atan2(Stacki[Nbval*Size+SP-1], Stacki[Nbval*Size+SP]);
 3111                        --SP; break;
 3112 
 3113           case cAtanh:
 3114             for(Nbval=0; Nbval<NbStack; Nbval++)
 3115             {
 3116               if(IsComplexType<Value_t>::result
 3117               ?  (Stacki[Nbval*Size+SP] == Value_t(-1) || Stacki[Nbval*Size+SP] == Value_t(1))
 3118               :  (Stacki[Nbval*Size+SP] <= Value_t(-1) || Stacki[Nbval*Size+SP] >= Value_t(1)))
 3119               { mData->mEvalErrorType=4; return Value_t(0); }
 3120                 Stacki[Nbval*Size+SP] = fp_atanh(Stacki[Nbval*Size+SP]);
 3121             }
 3122             break;
 3123 
 3124           case  cCbrt:
 3125             for(Nbval=0; Nbval<NbStack; Nbval++)
 3126                 Stacki[Nbval*Size+SP] = fp_cbrt(Stacki[Nbval*Size+SP]);
 3127             break;
 3128 
 3129           case  cCeil:
 3130             for(Nbval=0; Nbval<NbStack; Nbval++)
 3131                 Stacki[Nbval*Size+SP] = fp_ceil(Stacki[Nbval*Size+SP]);
 3132             break;
 3133 
 3134           case  cCosh:
 3135             for(Nbval=0; Nbval<NbStack; Nbval++)
 3136                 Stacki[Nbval*Size+SP] = fp_cosh(Stacki[Nbval*Size+SP]);
 3137             break;
 3138 
 3139           case   cCot:
 3140             for(Nbval=0; Nbval<NbStack; Nbval++)
 3141               {
 3142                   const Value_t t = fp_tan(Stacki[Nbval*Size+SP]);
 3143                   if(t == Value_t(0))
 3144                   { mData->mEvalErrorType=DIVISION_BY_ZERO; return Value_t(DIVISION_BY_ZERO); }
 3145                   Stacki[Nbval*Size+SP] = Value_t(1)/t;
 3146               }
 3147              break;
 3148 
 3149           case   cCsc:
 3150             for(Nbval=0; Nbval<NbStack; Nbval++)
 3151               {
 3152                   const Value_t s = fp_sin(Stacki[Nbval*Size+SP]);
 3153                   if(s == Value_t(0))
 3154                   { mData->mEvalErrorType=DIVISION_BY_ZERO; return Value_t(DIVISION_BY_ZERO); }
 3155                   Stacki[Nbval*Size+SP] = Value_t(1)/s;
 3156               }
 3157              break;
 3158 
 3159           case   cExp:
 3160             for(Nbval=0; Nbval<NbStack; Nbval++)
 3161                 Stacki[Nbval*Size+SP] = fp_exp(Stacki[Nbval*Size+SP]);
 3162             break;
 3163 
 3164           case   cExp2:
 3165             for(Nbval=0; Nbval<NbStack; Nbval++)
 3166                 Stacki[Nbval*Size+SP] = fp_exp2(Stacki[Nbval*Size+SP]);
 3167             break;
 3168 
 3169           case cFloor:
 3170             for(Nbval=0; Nbval<NbStack; Nbval++)
 3171                 Stacki[Nbval*Size+SP] = fp_floor(Stacki[Nbval*Size+SP]);
 3172             break;
 3173 
 3174           case cHypot:
 3175             for(Nbval=0; Nbval<NbStack; Nbval++)
 3176                 Stacki[Nbval*Size+SP-1] = fp_hypot(Stacki[Nbval*Size+SP-1], Stacki[Nbval*Size+SP]);
 3177             --SP;
 3178             break;
 3179 
 3180         case    cIf:
 3181         {
 3182             bool valid = fp_truth(Stacki[SP]);
 3183 
 3184             for(Nbval=1; Nbval<NbStack; Nbval++)
 3185                 if(valid != fp_truth(Stacki[Nbval*Size+SP]))
 3186                 {
 3187                     mData->mEvalErrorType=IF_FUNCT_ERROR; return Value_t(IF_FUNCT_ERROR);
 3188                 }
 3189             if(valid)
 3190                 IP += 2;
 3191             else
 3192             {
 3193                 const unsigned* buf = &byteCode[IP+1];
 3194                 IP = buf[0];
 3195                 DP = buf[1];
 3196             }
 3197         }
 3198             SP--;
 3199             break;
 3200       case   cInt:
 3201         for(Nbval=0; Nbval<NbStack; Nbval++)
 3202             Stacki[Nbval*Size+SP] = fp_int(Stacki[Nbval*Size+SP]);
 3203         break;
 3204 
 3205       case   cLog:
 3206         for(Nbval=0; Nbval<NbStack; Nbval++)
 3207         {
 3208           if(IsComplexType<Value_t>::result
 3209            ?   Stacki[Nbval*Size+SP] == Value_t(0)
 3210            :   !(Stacki[Nbval*Size+SP] > Value_t(0)))
 3211           { mData->mEvalErrorType=3; return Value_t(0); }
 3212               Stacki[Nbval*Size+SP] = fp_log(Stacki[Nbval*Size+SP]);
 3213         }
 3214         break;
 3215 
 3216           case cLog10:
 3217             for(Nbval=0; Nbval<NbStack; Nbval++)
 3218             {
 3219               if(IsComplexType<Value_t>::result
 3220                ?   Stacki[Nbval*Size+SP] == Value_t(0)
 3221                :   !(Stacki[Nbval*Size+SP] > Value_t(0)))
 3222               { mData->mEvalErrorType=3; return Value_t(0); }
 3223                   Stacki[Nbval*Size+SP] = fp_log10(Stacki[Nbval*Size+SP]);
 3224             }
 3225               break;
 3226 
 3227           case  cLog2:
 3228             for(Nbval=0; Nbval<NbStack; Nbval++)
 3229             {
 3230               if(IsComplexType<Value_t>::result
 3231                ?   Stacki[Nbval*Size+SP] == Value_t(0)
 3232                :   !(Stacki[Nbval*Size+SP] > Value_t(0)))
 3233               { mData->mEvalErrorType=3; return Value_t(0); }
 3234                   Stacki[Nbval*Size+SP] = fp_log2(Stacki[Nbval*Size+SP]);
 3235             }
 3236               break;
 3237 
 3238           case   cMax:
 3239             for(Nbval=0; Nbval<NbStack; Nbval++)
 3240                 Stacki[Nbval*Size+SP-1] = fp_max(Stacki[Nbval*Size+SP-1], Stacki[Nbval*Size+SP]);
 3241                        --SP; break;
 3242 
 3243           case   cMin:
 3244             for(Nbval=0; Nbval<NbStack; Nbval++)
 3245                 Stacki[Nbval*Size+SP-1] = fp_min(Stacki[Nbval*Size+SP-1], Stacki[Nbval*Size+SP]);
 3246                        --SP; break;
 3247 
 3248           case   cPow:
 3249               // x:Negative ^ y:NonInteger is failure,
 3250               // except when the reciprocal of y forms an integer
 3251               /*if(IsComplexType<Value_t>::result == false
 3252               && Stack[SP-1] < Value_t(0) &&
 3253                  !isInteger(Stack[SP]) &&
 3254                  !isInteger(1.0 / Stack[SP]))
 3255               { mEvalErrorType=3; return Value_t(0); }*/
 3256               // x:0 ^ y:negative is failure
 3257 
 3258             for(Nbval=0; Nbval<NbStack; Nbval++)
 3259             {
 3260               if(Stacki[Nbval*Size+SP-1] == Value_t(0) &&
 3261                  Stacki[Nbval*Size+SP] < Value_t(0))
 3262               { mData->mEvalErrorType=3; return Value_t(0); }
 3263                   Stacki[Nbval*Size+SP-1] = fp_pow(Stacki[Nbval*Size+SP-1], Stacki[Nbval*Size+SP]);
 3264             }
 3265               --SP; break;
 3266         case cPsh:
 3267           for(Nbval=0; Nbval<NbStack; Nbval++)
 3268           {
 3269               if((stt= NbStack*(Stacki[Nbval*Size+SP-1])+Nbval) >= VarSize || stt <0)
 3270               {
 3271                   mData->mEvalErrorType=VAR_OVERFLOW;
 3272                   return Value_t(VAR_OVERFLOW);
 3273               }
 3274               StackSave[stt] = Stacki[Nbval*Size+SP];
 3275               Stacki[Nbval*Size+SP-1] = 1.0;
 3276           }
 3277           --SP;
 3278           break;
 3279         case cCsd:
 3280           for(Nbval=0; Nbval<NbStack; Nbval++)
 3281           {
 3282               if((stt= NbStack*(Stacki[Nbval*Size+SP])+Nbval) >= VarSize || stt <0)
 3283               {
 3284                   mData->mEvalErrorType=VAR_OVERFLOW;
 3285                   return Value_t(VAR_OVERFLOW);
 3286               }
 3287               Stacki[Nbval*Size+SP] = StackSave[stt];
 3288           }
 3289           break;
 3290           case  cTrunc:
 3291             for(Nbval=0; Nbval<NbStack; Nbval++)
 3292                 Stacki[Nbval*Size+SP] = fp_trunc(Stacki[Nbval*Size+SP]);
 3293             break;
 3294 
 3295           case   cSec:
 3296             for(Nbval=0; Nbval<NbStack; Nbval++)
 3297               {
 3298                   const Value_t c = fp_cos(Stacki[Nbval*Size+SP]);
 3299                   if(c == Value_t(0))
 3300                   {
 3301                       mData->mEvalErrorType=DIVISION_BY_ZERO;
 3302                       return Value_t(DIVISION_BY_ZERO);
 3303                   }
 3304                   Stacki[Nbval*Size+SP] = Value_t(1)/c;
 3305               }
 3306               break;
 3307           case   cSin:
 3308             for(Nbval=0; Nbval<NbStack; Nbval++)
 3309                 Stacki[Nbval*Size+SP] = fp_sin(Stacki[Nbval*Size+SP]);
 3310             break;
 3311 
 3312           case  cSinh:
 3313             for(Nbval=0; Nbval<NbStack; Nbval++)
 3314                 Stacki[Nbval*Size+SP] = fp_sinh(Stacki[Nbval*Size+SP]);
 3315             break;
 3316 
 3317           case  cSqrt:
 3318             for(Nbval=0; Nbval<NbStack; Nbval++)
 3319             {
 3320               if(IsComplexType<Value_t>::result == false &&
 3321                  Stacki[Nbval*Size+SP] < Value_t(0))
 3322               {
 3323                   mData->mEvalErrorType=2;
 3324                   return Value_t(0);
 3325               }
 3326               Stacki[Nbval*Size+SP] = fp_sqrt(Stacki[Nbval*Size+SP]);
 3327             }
 3328             break;
 3329 
 3330           case   cTan:
 3331             for(Nbval=0; Nbval<NbStack; Nbval++)
 3332                 Stacki[Nbval*Size+SP] = fp_tan(Stacki[Nbval*Size+SP]);
 3333             break;
 3334 
 3335           case  cTanh:
 3336             for(Nbval=0; Nbval<NbStack; Nbval++)
 3337                 Stacki[Nbval*Size+SP] = fp_tanh(Stacki[Nbval*Size+SP]);
 3338             break;
 3339 
 3340 
 3341 // Misc:
 3342           case cImmed:
 3343             ++SP;
 3344             for(Nbval=0; Nbval<NbStack; Nbval++)
 3345                 Stacki[Nbval*Size+SP] = immed[DP];
 3346             DP++;
 3347             break;
 3348           case  cJump:
 3349               {
 3350                   const unsigned* buf = &byteCode[IP+1];
 3351                   IP = buf[0];
 3352                   DP = buf[1];
 3353                   break;
 3354               }
 3355 // Operators:
 3356           case   cNeg:
 3357             for(Nbval=0; Nbval<NbStack; Nbval++)
 3358                 Stacki[Nbval*Size+SP] = -Stacki[Nbval*Size+SP];
 3359             break;
 3360           case   cAdd:
 3361             for(Nbval=0; Nbval<NbStack; Nbval++)
 3362                 Stacki[Nbval*Size+SP-1] += Stacki[Nbval*Size+SP];
 3363             --SP;
 3364             break;
 3365           case   cSub:
 3366             for(Nbval=0; Nbval<NbStack; Nbval++)
 3367                 Stacki[Nbval*Size+SP-1] -= Stacki[Nbval*Size+SP];
 3368             --SP;
 3369             break;
 3370           case   cMul:
 3371             for(Nbval=0; Nbval<NbStack; Nbval++)
 3372                 Stacki[Nbval*Size+SP-1] *= Stacki[Nbval*Size+SP];
 3373             --SP;
 3374             break;
 3375           case   cDiv:
 3376             for(Nbval=0; Nbval<NbStack; Nbval++)
 3377             {
 3378               if(Stacki[Nbval*Size+SP] == Value_t(0))
 3379               {
 3380                   mData->mEvalErrorType=DIVISION_BY_ZERO;
 3381                   return Value_t(DIVISION_BY_ZERO);
 3382               }
 3383               Stacki[Nbval*Size+SP-1] /= Stacki[Nbval*Size+SP];
 3384             }
 3385             --SP; break;
 3386           case   cMod:
 3387             for(Nbval=0; Nbval<NbStack; Nbval++)
 3388             {
 3389               if(Stacki[Nbval*Size+SP] == Value_t(0))
 3390               {
 3391                   mData->mEvalErrorType=DIVISION_BY_ZERO;
 3392                   return Value_t(DIVISION_BY_ZERO);
 3393               }
 3394               Stacki[Nbval*Size+SP-1] = fp_mod(Stacki[Nbval*Size+SP-1], Stacki[Nbval*Size+SP]);
 3395             }
 3396             --SP;
 3397             break;
 3398           case cEqual:
 3399             for(Nbval=0; Nbval<NbStack; Nbval++)
 3400               Stacki[Nbval*Size+SP-1] = fp_equal(Stacki[Nbval*Size+SP-1], Stacki[Nbval*Size+SP]);
 3401             --SP;
 3402             break;
 3403           case cNEqual:
 3404             for(Nbval=0; Nbval<NbStack; Nbval++)
 3405               Stacki[Nbval*Size+SP-1] = fp_nequal(Stacki[Nbval*Size+SP-1], Stacki[Nbval*Size+SP]);
 3406             --SP;
 3407             break;
 3408           case  cLess:
 3409             for(Nbval=0; Nbval<NbStack; Nbval++)
 3410               Stacki[Nbval*Size+SP-1] = fp_less(Stacki[Nbval*Size+SP-1], Stacki[Nbval*Size+SP]);
 3411             --SP;
 3412             break;
 3413 
 3414           case  cLessOrEq:
 3415             for(Nbval=0; Nbval<NbStack; Nbval++)
 3416               Stacki[Nbval*Size+SP-1] = fp_lessOrEq(Stacki[Nbval*Size+SP-1], Stacki[Nbval*Size+SP]);
 3417             --SP;
 3418             break;
 3419 
 3420           case cGreater:
 3421             for(Nbval=0; Nbval<NbStack; Nbval++)
 3422               Stacki[Nbval*Size+SP-1] = fp_less(Stacki[Nbval*Size+SP], Stacki[Nbval*Size+SP-1]);
 3423             --SP;
 3424             break;
 3425 
 3426           case cGreaterOrEq:
 3427             for(Nbval=0; Nbval<NbStack; Nbval++)
 3428               Stacki[Nbval*Size+SP-1] = fp_lessOrEq(Stacki[Nbval*Size+SP], Stacki[Nbval*Size+SP-1]);
 3429             --SP;
 3430             break;
 3431 
 3432           case   cNot:
 3433             for(Nbval=0; Nbval<NbStack; Nbval++)
 3434                 Stacki[Nbval*Size+SP] = fp_not(Stacki[Nbval*Size+SP]);
 3435             break;
 3436 
 3437           case cNotNot:
 3438             for(Nbval=0; Nbval<NbStack; Nbval++)
 3439                 Stacki[Nbval*Size+SP] = fp_notNot(Stacki[Nbval*Size+SP]);
 3440             break;
 3441 
 3442           case   cAnd:
 3443             for(Nbval=0; Nbval<NbStack; Nbval++)
 3444               Stacki[Nbval*Size+SP-1] = fp_and(Stacki[Nbval*Size+SP-1], Stacki[Nbval*Size+SP]);
 3445             --SP;
 3446             break;
 3447 
 3448           case    cOr:
 3449             for(Nbval=0; Nbval<NbStack; Nbval++)
 3450               Stacki[Nbval*Size+SP-1] = fp_or(Stacki[Nbval*Size+SP-1], Stacki[Nbval*Size+SP]);
 3451             --SP;
 3452             break;
 3453 
 3454 // Degrees-radians conversion:
 3455           case   cDeg:
 3456             for(Nbval=0; Nbval<NbStack; Nbval++)
 3457                 Stacki[Nbval*Size+SP] = RadiansToDegrees(Stacki[Nbval*Size+SP]);
 3458             break;
 3459           case   cRad:
 3460             for(Nbval=0; Nbval<NbStack; Nbval++)
 3461                 Stacki[Nbval*Size+SP] = DegreesToRadians(Stacki[Nbval*Size+SP]);
 3462             break;
 3463 
 3464 // User-defined function calls:
 3465           case cFCall:
 3466               {
 3467                   const unsigned index = byteCode[++IP];
 3468                   const unsigned params = mData->mFuncPtrs[index].mParams;
 3469                   for(Nbval=0; Nbval<NbStack; Nbval++)
 3470                   {
 3471                   const Value_t retVal =
 3472                       mData->mFuncPtrs[index].mRawFuncPtr ?
 3473                       mData->mFuncPtrs[index].mRawFuncPtr(&(Stacki[Nbval*Size+SP-params+1])) :
 3474                       mData->mFuncPtrs[index].mFuncWrapperPtr->callFunction
 3475                       (&Stacki[Nbval*Size+SP-params+1]);
 3476                   Stacki[Nbval*Size+SP-(int(params)-1)] = retVal;
 3477                   }
 3478                   SP -= int(params)-1;
 3479                   break;
 3480               }
 3481 
 3482           case   cFetch:
 3483               {
 3484                   unsigned stackOffs = byteCode[++IP];
 3485                   for(Nbval=0; Nbval<NbStack; Nbval++)
 3486                     Stacki[Nbval*Size+SP+1] = Stacki[Nbval*Size+stackOffs];
 3487                   ++SP;
 3488                   break;
 3489               }
 3490 
 3491 #ifdef FP_SUPPORT_OPTIMIZER
 3492           case   cPopNMov:
 3493               {
 3494                   unsigned stackOffs_target = byteCode[++IP];
 3495                   unsigned stackOffs_source = byteCode[++IP];
 3496                   for(Nbval=0; Nbval<NbStack; Nbval++)
 3497                       Stacki[Nbval*Size+stackOffs_target] = Stacki[Nbval*Size+stackOffs_source];
 3498                   SP = stackOffs_target;
 3499                   break;
 3500               }
 3501 
 3502           case  cLog2by:
 3503             for(Nbval=0; Nbval<NbStack; Nbval++)
 3504             {
 3505               if(IsComplexType<Value_t>::result
 3506                ?   Stacki[Nbval*Size+SP-1] == Value_t(0)
 3507                :   !(Stacki[Nbval*Size+SP-1] > Value_t(0)))
 3508               { mData->mEvalErrorType=3; return Value_t(0); }
 3509                   Stacki[Nbval*Size+SP-1] = fp_log2(Stacki[Nbval*Size+SP-1]) * Stacki[Nbval*Size+SP];
 3510             }
 3511               --SP;
 3512               break;
 3513 
 3514           case cNop: break;
 3515 #endif // FP_SUPPORT_OPTIMIZER
 3516 
 3517           case cSinCos:
 3518               for(Nbval=0; Nbval<NbStack; Nbval++)
 3519                   fp_sinCos(Stacki[Nbval*Size+SP], Stacki[Nbval*Size+SP+1], Stacki[Nbval*Size+SP]);
 3520               ++SP;
 3521               break;
 3522           case cSinhCosh:
 3523               for(Nbval=0; Nbval<NbStack; Nbval++)
 3524                   fp_sinhCosh(Stacki[Nbval*Size+SP], Stacki[Nbval*Size+SP+1], Stacki[Nbval*Size+SP]);
 3525               ++SP;
 3526               break;
 3527 
 3528           case cAbsNot:
 3529               for(Nbval=0; Nbval<NbStack; Nbval++)
 3530                   Stacki[Nbval*Size+SP] = fp_absNot(Stacki[Nbval*Size+SP]);
 3531             break;
 3532           case cAbsNotNot:
 3533               for(Nbval=0; Nbval<NbStack; Nbval++)
 3534                   Stacki[Nbval*Size+SP] = fp_absNotNot(Stacki[Nbval*Size+SP]);
 3535             break;
 3536           case cAbsAnd:
 3537               for(Nbval=0; Nbval<NbStack; Nbval++)
 3538                   Stacki[Nbval*Size+SP-1] = fp_absAnd(Stacki[Nbval*Size+SP-1], Stacki[Nbval*Size+SP]);
 3539               --SP; break;
 3540           case cAbsOr:
 3541               for(Nbval=0; Nbval<NbStack; Nbval++)
 3542                   Stacki[Nbval*Size+SP-1] = fp_absOr(Stacki[Nbval*Size+SP-1], Stacki[Nbval*Size+SP]);
 3543               --SP; break;
 3544           case cAbsIf:
 3545             {
 3546                 bool valid = fp_absTruth(Stacki[SP]);
 3547 
 3548                 for(Nbval=1; Nbval<NbStack; Nbval++)
 3549                     if(valid != fp_absTruth(Stacki[Nbval*Size+SP]))
 3550                         { mData->mEvalErrorType=IF_FUNCT_ERROR; return Value_t(IF_FUNCT_ERROR); }
 3551 
 3552                 if(valid)
 3553                     IP += 2;
 3554                 else
 3555                 {
 3556                     const unsigned* buf = &byteCode[IP+1];
 3557                     IP = buf[0];
 3558                     DP = buf[1];
 3559                 }
 3560             }
 3561             SP--;
 3562             break;
 3563           case   cDup:
 3564             for(Nbval=0; Nbval<NbStack; Nbval++)
 3565                 Stacki[Nbval*Size+SP+1] = Stacki[Nbval*Size+SP];
 3566             ++SP; break;
 3567 
 3568           case   cInv:
 3569             for(Nbval=0; Nbval<NbStack; Nbval++)
 3570             {
 3571               if(Stacki[Nbval*Size+SP] == Value_t(0))
 3572               { mData->mEvalErrorType=DIVISION_BY_ZERO; return Value_t(DIVISION_BY_ZERO); }
 3573               Stacki[Nbval*Size+SP] = Value_t(1)/Stacki[Nbval*Size+SP];
 3574             }
 3575               break;
 3576 
 3577           case   cSqr:
 3578               for(Nbval=0; Nbval<NbStack; Nbval++)
 3579                       Stacki[Nbval*Size+SP] = Stacki[Nbval*Size+SP]*Stacki[Nbval*Size+SP];
 3580               break;
 3581 
 3582           case   cRDiv:
 3583             for(Nbval=0; Nbval<NbStack; Nbval++)
 3584             {
 3585               if(Stacki[Nbval*Size+SP-1] == Value_t(0))
 3586               { mData->mEvalErrorType=DIVISION_BY_ZERO; return Value_t(DIVISION_BY_ZERO); }
 3587               Stacki[Nbval*Size+SP-1] = Stacki[Nbval*Size+SP] / Stacki[Nbval*Size+SP-1];
 3588             }
 3589             --SP; break;
 3590 
 3591           case   cRSub:
 3592             for(Nbval=0; Nbval<NbStack; Nbval++)
 3593                 Stacki[Nbval*Size+SP-1] = Stacki[Nbval*Size+SP] - Stacki[Nbval*Size+SP-1];
 3594             --SP; break;
 3595 
 3596           case   cRSqrt:
 3597               for(Nbval=0; Nbval<NbStack; Nbval++)
 3598               {
 3599                 if(Stacki[Nbval*Size+SP] == Value_t(0))
 3600                 { mData->mEvalErrorType=DIVISION_BY_ZERO; return Value_t(DIVISION_BY_ZERO); }
 3601                 Stacki[Nbval*Size+SP] = Value_t(1) / fp_sqrt(Stacki[Nbval*Size+SP]);
 3602               }
 3603                 break;
 3604 
 3605 #ifdef FP_SUPPORT_COMPLEX_NUMBERS
 3606           case   cReal:for(Nbval=0; Nbval<NbStack; Nbval++)
 3607                 Stacki[Nbval*Size+SP] = fp_real(Stacki[Nbval*Size+SP]); break;
 3608           case   cImag:for(Nbval=0; Nbval<NbStack; Nbval++)
 3609                 Stacki[Nbval*Size+SP] = fp_imag(Stacki[Nbval*Size+SP]); break;
 3610           case   cArg: for(Nbval=0; Nbval<NbStack; Nbval++)
 3611                 Stacki[Nbval*Size+SP] = fp_arg(Stacki[Nbval*Size+SP]); break;
 3612           case   cConj:for(Nbval=0; Nbval<NbStack; Nbval++)
 3613                 Stacki[Nbval*Size+SP] = fp_conj(Stacki[Nbval*Size+SP]); break;
 3614           case   cPolar:for(Nbval=0; Nbval<NbStack; Nbval++)
 3615                 Stacki[Nbval*Size+SP-1] = fp_polar(Stacki[Nbval*Size+SP-1], Stacki[Nbval*Size+SP]);
 3616               --SP; break;
 3617 #endif
 3618 
 3619         case cPCall:
 3620             {
 3621                 unsigned index = byteCode[++IP];
 3622                 unsigned params = mData->mFuncParsers[index].mParams;
 3623                 double res[NbStack];
 3624                 double rest=mData->mFuncParsers[index].mParserPtr->Eval2
 3625                 (&(Stacki[SP-params+1]), Size, res, NbStack);
 3626                 if (int(rest) == IF_FUNCT_ERROR)
 3627                 {
 3628                     mData->mEvalErrorType = IF_FUNCT_ERROR;
 3629                     return IF_FUNCT_ERROR;
 3630                 }
 3631                 if (int(rest) == VAR_OVERFLOW)
 3632                 {
 3633                     mData->mEvalErrorType = VAR_OVERFLOW;
 3634                     return VAR_OVERFLOW;
 3635                 }
 3636                 for(Nbval=0; Nbval<NbStack; Nbval++)
 3637                 {
 3638                   Stacki[Nbval*Size+SP - (int(params)-1)] = res[Nbval];
 3639                 }
 3640                 SP -= int(params)-1;
 3641                 break;
 3642             }
 3643 // Variables:
 3644           default:
 3645             ++SP;
 3646             for(Nbval=0; Nbval<NbStack; Nbval++)
 3647               Stacki[Nbval*Size+SP] = Vars[byteCode[IP]+ NbVar*Nbval-VarBegin];
 3648         }
 3649     }
 3650 
 3651     mData->mEvalErrorType=EVAL_NO_ERROR;
 3652     for(Nbval=0; Nbval<NbStack; Nbval++)
 3653         results[Nbval] = Stacki[Nbval*Size+SP];
 3654     return Value_t(EVAL_NO_ERROR);
 3655 }
 3656 
 3657 
 3658 //===========================================================================
 3659 // Variable deduction
 3660 //===========================================================================
 3661 namespace
 3662 {
 3663     template<typename Value_t>
 3664     int deduceVariables(FunctionParserBase<Value_t>& fParser,
 3665                         const char* funcStr,
 3666                         std::string& destVarString,
 3667                         int* amountOfVariablesFound,
 3668                         std::vector<std::string>* destVarNames,
 3669                         bool useDegrees)
 3670     {
 3671         typedef std::set<std::string> StrSet;
 3672         StrSet varNames;
 3673 
 3674         int oldIndex = -1;
 3675 
 3676         while(true)
 3677         {
 3678             destVarString.clear();
 3679             for(StrSet::iterator iter = varNames.begin();
 3680                 iter != varNames.end();
 3681                 ++iter)
 3682             {
 3683                 if(iter != varNames.begin()) destVarString += ",";
 3684                 destVarString += *iter;
 3685             }
 3686 
 3687             const int index =
 3688                 fParser.Parse(funcStr, destVarString, useDegrees);
 3689             if(index < 0) break;
 3690             if(index == oldIndex) return index;
 3691 
 3692             unsigned nameLength = readIdentifier<Value_t>(funcStr + index);
 3693             if(nameLength & 0x80000000U) return index;
 3694             if(nameLength == 0) return index;
 3695 
 3696             varNames.insert(std::string(funcStr + index, nameLength));
 3697             oldIndex = index;
 3698         }
 3699 
 3700         if(amountOfVariablesFound)
 3701             *amountOfVariablesFound = int(varNames.size());
 3702 
 3703         if(destVarNames)
 3704             destVarNames->assign(varNames.begin(), varNames.end());
 3705 
 3706         return -1;
 3707     }
 3708 }
 3709 
 3710 template<typename Value_t>
 3711 int FunctionParserBase<Value_t>::ParseAndDeduceVariables
 3712 (const std::string& function,
 3713  int* amountOfVariablesFound,
 3714  bool useDegrees)
 3715 {
 3716     std::string varString;
 3717     return deduceVariables(*this, function.c_str(), varString,
 3718                            amountOfVariablesFound, 0, useDegrees);
 3719 }
 3720 
 3721 template<typename Value_t>
 3722 int FunctionParserBase<Value_t>::ParseAndDeduceVariables
 3723 (const std::string& function,
 3724  std::string& resultVarString,
 3725  int* amountOfVariablesFound,
 3726  bool useDegrees)
 3727 {
 3728     std::string varString;
 3729     const int index =
 3730         deduceVariables(*this, function.c_str(), varString,
 3731                         amountOfVariablesFound, 0, useDegrees);
 3732     if(index < 0) resultVarString = varString;
 3733     return index;
 3734 }
 3735 
 3736 template<typename Value_t>
 3737 int FunctionParserBase<Value_t>::ParseAndDeduceVariables
 3738 (const std::string& function,
 3739  std::vector<std::string>& resultVars,
 3740  bool useDegrees)
 3741 {
 3742     std::string varString;
 3743     std::vector<std::string> vars;
 3744     const int index =
 3745         deduceVariables(*this, function.c_str(), varString,
 3746                         0, &vars, useDegrees);
 3747     if(index < 0) resultVars.swap(vars);
 3748     return index;
 3749 }
 3750 
 3751 
 3752 #ifdef FUNCTIONPARSER_SUPPORT_DEBUGGING
 3753 //===========================================================================
 3754 // Bytecode injection
 3755 //===========================================================================
 3756 template<typename Value_t>
 3757 void FunctionParserBase<Value_t>::InjectRawByteCode
 3758 (const unsigned* bytecode, unsigned bytecodeAmount,
 3759  const Value_t* immed, unsigned immedAmount, unsigned stackSize)
 3760 {
 3761     CopyOnWrite();
 3762 
 3763     mData->mByteCode.assign(bytecode, bytecode + bytecodeAmount);
 3764     mData->mImmed.assign(immed, immed + immedAmount);
 3765     mData->mStackSize = stackSize;
 3766 
 3767 #ifndef FP_USE_THREAD_SAFE_EVAL
 3768     mData->mStack.resize(stackSize);
 3769     //mData->mStacki.resize(64*mData->mStackSize);
 3770 #endif
 3771 }
 3772 
 3773 //===========================================================================
 3774 // Debug output
 3775 //===========================================================================
 3776 #include <iomanip>
 3777 #include <sstream>
 3778 namespace
 3779 {
 3780     inline void printHex(std::ostream& dest, unsigned n)
 3781     {
 3782         std::ios::fmtflags flags = dest.flags();
 3783         dest.width(4); dest.fill('0'); std::hex(dest); //uppercase(dest);
 3784         dest << n;
 3785         dest.flags(flags);
 3786     }
 3787 
 3788     void padLine(std::ostringstream& dest, unsigned destLength)
 3789     {
 3790         for(std::size_t currentLength = dest.str().length();
 3791             currentLength < destLength;
 3792             ++currentLength)
 3793         {
 3794             dest << ' ';
 3795         }
 3796     }
 3797 
 3798     const struct PowiMuliType
 3799     {
 3800         unsigned opcode_square;
 3801         unsigned opcode_cumulate;
 3802         unsigned opcode_invert;
 3803         unsigned opcode_half;
 3804         unsigned opcode_invhalf;
 3805     } iseq_powi = {cSqr,cMul,cInv,cSqrt,cRSqrt},
 3806       iseq_muli = {~unsigned(0), cAdd,cNeg, ~unsigned(0),~unsigned(0) };
 3807 
 3808     template<typename Value_t>
 3809     Value_t ParsePowiMuli(
 3810         const PowiMuliType& opcodes,
 3811         const std::vector<unsigned>& ByteCode, unsigned& IP,
 3812         unsigned limit,
 3813         std::size_t factor_stack_base,
 3814         std::vector<Value_t>& stack,
 3815         bool IgnoreExcess)
 3816     {
 3817         Value_t result = Value_t(1);
 3818         while(IP < limit)
 3819         {
 3820             if(ByteCode[IP] == opcodes.opcode_square)
 3821             {
 3822                 if(!isInteger(result)) break;
 3823                 result *= Value_t(2);
 3824                 ++IP;
 3825                 continue;
 3826             }
 3827             if(ByteCode[IP] == opcodes.opcode_invert)
 3828             {
 3829                 if(result < Value_t(0)) break;
 3830                 result = -result;
 3831                 ++IP;
 3832                 continue;
 3833             }
 3834             if(ByteCode[IP] == opcodes.opcode_half)
 3835             {
 3836                 if(result > Value_t(0) && isEvenInteger(result))
 3837                     break;
 3838                 if(isInteger(result * Value_t(0.5))) break;
 3839                 result *= Value_t(0.5);
 3840                 ++IP;
 3841                 continue;
 3842             }
 3843             if(ByteCode[IP] == opcodes.opcode_invhalf)
 3844             {
 3845                 if(result > Value_t(0) && isEvenInteger(result))
 3846                     break;
 3847                 if(isInteger(result * Value_t(-0.5))) break;
 3848                 result *= Value_t(-0.5);
 3849                 ++IP;
 3850                 continue;
 3851             }
 3852 
 3853             unsigned dup_fetch_pos = IP;
 3854             Value_t lhs = Value_t(1);
 3855 
 3856             if(ByteCode[IP] == cFetch)
 3857             {
 3858                 unsigned index = ByteCode[++IP];
 3859                 if(index < factor_stack_base
 3860                 || std::size_t(index-factor_stack_base) >= stack.size())
 3861                 {
 3862                     // It wasn't a powi-fetch after all
 3863                     IP = dup_fetch_pos;
 3864                     break;
 3865                 }
 3866                 lhs = stack[index - factor_stack_base];
 3867                 // Note: ^This assumes that cFetch of recentmost
 3868                 //        is always converted into cDup.
 3869                 goto dup_or_fetch;
 3870             }
 3871 
 3872             if(ByteCode[IP] == cDup)
 3873             {
 3874                 lhs = result;
 3875                 goto dup_or_fetch;
 3876 
 3877             dup_or_fetch:
 3878                 stack.push_back(result);
 3879                 ++IP;
 3880                 Value_t subexponent = ParsePowiMuli
 3881                     (opcodes,
 3882                      ByteCode, IP, limit,
 3883                      factor_stack_base, stack,
 3884                      IgnoreExcess);
 3885                 if(IP >= limit && IgnoreExcess)
 3886                     return lhs*subexponent;
 3887                 if(IP >= limit || ByteCode[IP] != opcodes.opcode_cumulate)
 3888                 {
 3889                     // It wasn't a powi-dup after all
 3890                     IP = dup_fetch_pos;
 3891                     break;
 3892                 }
 3893                 ++IP; // skip opcode_cumulate
 3894                 stack.pop_back();
 3895                 result += lhs*subexponent;
 3896                 continue;
 3897             }
 3898             break;
 3899         }
 3900         return result;
 3901     }
 3902 
 3903     template<typename Value_t>
 3904     Value_t ParsePowiSequence(const std::vector<unsigned>& ByteCode,
 3905                               unsigned& IP, unsigned limit,
 3906                               std::size_t factor_stack_base,
 3907                               bool IgnoreExcess = false)
 3908     {
 3909         std::vector<Value_t> stack;
 3910         stack.push_back(Value_t(1));
 3911         return ParsePowiMuli(iseq_powi, ByteCode, IP, limit,
 3912                              factor_stack_base, stack,
 3913                              IgnoreExcess);
 3914     }
 3915 
 3916     template<typename Value_t>
 3917     Value_t ParseMuliSequence(const std::vector<unsigned>& ByteCode,
 3918                               unsigned& IP, unsigned limit,
 3919                               std::size_t factor_stack_base,
 3920                               bool IgnoreExcess = false)
 3921     {
 3922         std::vector<Value_t> stack;
 3923         stack.push_back(Value_t(1));
 3924         return ParsePowiMuli(iseq_muli, ByteCode, IP, limit,
 3925                              factor_stack_base, stack,
 3926                              IgnoreExcess);
 3927     }
 3928 
 3929     struct IfInfo
 3930     {
 3931         std::pair<int,std::string> condition;
 3932         std::pair<int,std::string> thenbranch;
 3933         unsigned endif_location;
 3934 
 3935         IfInfo() : condition(), thenbranch(), endif_location() { }
 3936     };
 3937 }
 3938 
 3939 template<typename Value_t>
 3940 void FunctionParserBase<Value_t>::PrintByteCode(std::ostream& dest,
 3941                                                 bool showExpression) const
 3942 {
 3943     dest << "Size of stack: " << mData->mStackSize << "\n";
 3944 
 3945     std::ostringstream outputBuffer;
 3946     std::ostream& output = (showExpression ? outputBuffer : dest);
 3947 
 3948     const std::vector<unsigned>& ByteCode = mData->mByteCode;
 3949     const std::vector<Value_t>& Immed = mData->mImmed;
 3950 
 3951     std::vector<std::pair<int,std::string> > stack;
 3952     std::vector<IfInfo> if_stack;
 3953 
 3954     for(unsigned IP = 0, DP = 0; IP <= ByteCode.size(); ++IP)
 3955     {
 3956     after_powi_or_muli:;
 3957         std::string n;
 3958         bool out_params = false;
 3959         unsigned params = 2, produces = 1, opcode = 0;
 3960 
 3961         if(showExpression && !if_stack.empty() &&
 3962           (   // Normal If termination rule:
 3963               if_stack.back().endif_location == IP
 3964               // This rule matches when cJumps are threaded:
 3965            || (IP < ByteCode.size() && ByteCode[IP] == cJump
 3966                && !if_stack.back().thenbranch.second.empty())
 3967           ))
 3968         {
 3969             printHex(output, IP);
 3970             if(if_stack.back().endif_location == IP)
 3971                 output << ": ----- (phi)";
 3972             else
 3973                 output << ": ----- (phi+)";
 3974 
 3975             stack.resize(stack.size()+2);
 3976             std::swap(stack[stack.size()-3], stack[stack.size()-1]);
 3977             std::swap(if_stack.back().condition,  stack[stack.size()-3]);
 3978             std::swap(if_stack.back().thenbranch, stack[stack.size()-2]);
 3979             opcode = cIf;
 3980             params = 3;
 3981             --IP;
 3982             if_stack.pop_back();
 3983         }
 3984         else
 3985         {
 3986             if(IP >= ByteCode.size()) break;
 3987             opcode = ByteCode[IP];
 3988 
 3989             if(showExpression && (
 3990                 opcode == cSqr || opcode == cDup
 3991              || opcode == cInv
 3992              || opcode == cSqrt || opcode == cRSqrt
 3993              || opcode == cFetch
 3994             ))
 3995             {
 3996                 unsigned changed_ip = IP;
 3997                 Value_t exponent =
 3998                     ParsePowiSequence<Value_t>
 3999                     (ByteCode, changed_ip,
 4000                      if_stack.empty()
 4001                      ? (unsigned)ByteCode.size()
 4002                      : if_stack.back().endif_location,
 4003                      stack.size()-1);
 4004                 std::string        operation_prefix;
 4005                 std::ostringstream operation_value;
 4006                 int prio = 0;
 4007                 if(exponent == Value_t(1.0))
 4008                 {
 4009                     if(opcode != cDup) goto not_powi_or_muli;
 4010                     Value_t factor =
 4011                         ParseMuliSequence<Value_t>
 4012                         (ByteCode, changed_ip,
 4013                          if_stack.empty()
 4014                          ? (unsigned)ByteCode.size()
 4015                          : if_stack.back().endif_location,
 4016                          stack.size()-1);
 4017                     if(factor == Value_t(1) || factor == Value_t(-1))
 4018                         goto not_powi_or_muli;
 4019                     operation_prefix = "*";
 4020                     operation_value << factor;
 4021                     prio = 3;
 4022                 }
 4023                 else
 4024                 {
 4025                     prio = 2;
 4026                     operation_prefix = "^";
 4027                     operation_value << exponent;
 4028                 }
 4029 
 4030                 //unsigned explanation_before = changed_ip-2;
 4031                 unsigned explanation_before = changed_ip-1;
 4032 
 4033                 const char* explanation_prefix = "_";
 4034                 for(const unsigned first_ip = IP; IP < changed_ip; ++IP)
 4035                 {
 4036                     printHex(output, IP);
 4037                     output << ": ";
 4038 
 4039                     const char* sep = "|";
 4040                     if(first_ip+1 == changed_ip)
 4041                     { sep = "="; explanation_prefix = " "; }
 4042                     else if(IP   == first_ip) sep = "\\";
 4043                     else if(IP+1 == changed_ip) sep = "/";
 4044                     else explanation_prefix = "=";
 4045 
 4046                     switch(ByteCode[IP])
 4047                     {
 4048                         case cInv: output << "inv"; break;
 4049                         case cNeg: output << "neg"; break;
 4050                         case cDup: output << "dup"; break;
 4051                         case cSqr: output << "sqr"; break;
 4052                         case cMul: output << "mul"; break;
 4053                         case cAdd: output << "add"; break;
 4054                         case cCbrt: output << "cbrt"; break;
 4055                         case cSqrt: output << "sqrt"; break;
 4056                         case cRSqrt: output << "rsqrt"; break;
 4057                         case cFetch:
 4058                         {
 4059                             unsigned index = ByteCode[++IP];
 4060                             output << "cFetch(" << index << ")";
 4061                             break;
 4062                         }
 4063                         default: break;
 4064                     }
 4065                     padLine(outputBuffer, 20);
 4066                     output << sep;
 4067                     if(IP >= explanation_before)
 4068                     {
 4069                         explanation_before = (unsigned)ByteCode.size();
 4070                         output << explanation_prefix
 4071                                << '[' << (stack.size()-1) << ']';
 4072                         std::string last = stack.back().second;
 4073                         if(stack.back().first >= prio)
 4074                             last = "(" + last + ")";
 4075                         output << last;
 4076                         output << operation_prefix;
 4077                         output << operation_value.str();
 4078                     }
 4079                     else
 4080                     {
 4081                         unsigned p = first_ip;
 4082                         Value_t exp = operation_prefix=="^" ?
 4083                             ParsePowiSequence<Value_t>
 4084                             (ByteCode, p, IP+1, stack.size()-1, true) :
 4085                             ParseMuliSequence<Value_t>
 4086                             (ByteCode, p, IP+1, stack.size()-1, true);
 4087                         std::string last = stack.back().second;
 4088                         if(stack.back().first >= prio)
 4089                             last = "(" + last + ")";
 4090                         output << " ..." << last;
 4091                         output << operation_prefix;
 4092                         output << exp;
 4093                     }
 4094                     dest << outputBuffer.str() << std::endl;
 4095                     outputBuffer.str("");
 4096                 }
 4097 
 4098                 std::string& last = stack.back().second;
 4099                 if(stack.back().first >= prio)
 4100                     last = "(" + last + ")";
 4101                 last += operation_prefix;
 4102                 last += operation_value.str();
 4103                 stack.back().first = prio;
 4104 
 4105                 goto after_powi_or_muli;
 4106             }
 4107         not_powi_or_muli:;
 4108             printHex(output, IP);
 4109             output << ": ";
 4110 
 4111             switch(opcode)
 4112             {
 4113               case cIf:
 4114               {
 4115                   unsigned label = ByteCode[IP+1]+1;
 4116                   output << "jz ";
 4117                   printHex(output, label);
 4118                   params = 1;
 4119                   produces = 0;
 4120                   IP += 2;
 4121 
 4122                   if_stack.resize(if_stack.size() + 1);
 4123                   std::swap( if_stack.back().condition, stack.back() );
 4124                   if_stack.back().endif_location = (unsigned) ByteCode.size();
 4125                   stack.pop_back();
 4126                   break;
 4127               }
 4128               case cAbsIf:
 4129               {
 4130                   unsigned dp    = ByteCode[IP+2];
 4131                   unsigned label = ByteCode[IP+1]+1;
 4132                   output << "jz_abs " << dp << ",";
 4133                   printHex(output, label);
 4134                   params = 1;
 4135                   produces = 0;
 4136                   IP += 2;
 4137 
 4138                   if_stack.resize(if_stack.size() + 1);
 4139                   std::swap( if_stack.back().condition, stack.back() );
 4140                   if_stack.back().endif_location = (unsigned) ByteCode.size();
 4141                   stack.pop_back();
 4142                   break;
 4143               }
 4144 
 4145               case cJump:
 4146               {
 4147                   unsigned dp    = ByteCode[IP+2];
 4148                   unsigned label = ByteCode[IP+1]+1;
 4149 
 4150                   if(!if_stack.empty() && !stack.empty())
 4151                   {
 4152                       std::swap(if_stack.back().thenbranch, stack.back());
 4153                       if_stack.back().endif_location = label;
 4154                       stack.pop_back();
 4155                   }
 4156 
 4157                   output << "jump " << dp << ",";
 4158                   printHex(output, label);
 4159                   params = 0;
 4160                   produces = 0;
 4161                   IP += 2;
 4162                   break;
 4163               }
 4164               case cImmed:
 4165               {
 4166                   if(showExpression)
 4167                   {
 4168                       std::ostringstream buf;
 4169                       buf.precision(8);
 4170                       buf << Immed[DP];
 4171                       stack.push_back( std::make_pair(0, buf.str()) );
 4172                   }
 4173                   output.precision(8);
 4174                   output << "push " << Immed[DP];
 4175                   ++DP;
 4176                   produces = 0;
 4177                   break;
 4178               }
 4179 
 4180               case cFCall:
 4181                   {
 4182                       const unsigned index = ByteCode[++IP];
 4183                       params = mData->mFuncPtrs[index].mParams;
 4184                       static std::string name;
 4185                       name = "f:" + findName(mData->mNamePtrs, index,
 4186                                              NameData<Value_t>::FUNC_PTR);
 4187                       n = name.c_str();
 4188                       out_params = true;
 4189                       break;
 4190                   }
 4191 
 4192               case cPCall:
 4193                   {
 4194                       const unsigned index = ByteCode[++IP];
 4195                       params = mData->mFuncParsers[index].mParams;
 4196                       static std::string name;
 4197                       name = "p:" + findName(mData->mNamePtrs, index,
 4198                                              NameData<Value_t>::PARSER_PTR);
 4199                       n = name.c_str();
 4200                       out_params = true;
 4201                       break;
 4202                   }
 4203 
 4204               default:
 4205                   if(IsVarOpcode(opcode))
 4206                   {
 4207                       if(showExpression)
 4208                       {
 4209                           stack.push_back(std::make_pair(0,
 4210                               (findName(mData->mNamePtrs, opcode,
 4211                                         NameData<Value_t>::VARIABLE))));
 4212                       }
 4213                       output << "push Var" << opcode-VarBegin;
 4214                       produces = 0;
 4215                   }
 4216                   else
 4217                   {
 4218                       switch(OPCODE(opcode))
 4219                       {
 4220                         case cNeg: n = "neg"; params = 1; break;
 4221                         case cAdd: n = "add"; break;
 4222                         case cSub: n = "sub"; break;
 4223                         case cMul: n = "mul"; break;
 4224                         case cDiv: n = "div"; break;
 4225                         case cMod: n = "mod"; break;
 4226                         case cPow: n = "pow"; break;
 4227                         case cEqual: n = "eq"; break;
 4228                         case cNEqual: n = "neq"; break;
 4229                         case cLess: n = "lt"; break;
 4230                         case cLessOrEq: n = "le"; break;
 4231                         case cGreater: n = "gt"; break;
 4232                         case cGreaterOrEq: n = "ge"; break;
 4233                         case cAnd: n = "and"; break;
 4234                         case cOr: n = "or"; break;
 4235                         case cNot: n = "not"; params = 1; break;
 4236                         case cNotNot: n = "notnot"; params = 1; break;
 4237                         case cDeg: n = "deg"; params = 1; break;
 4238                         case cRad: n = "rad"; params = 1; break;
 4239 
 4240                         case cFetch:
 4241                         {
 4242                             unsigned index = ByteCode[++IP];
 4243                             if(showExpression && index < stack.size())
 4244                                 stack.push_back(stack[index]);
 4245                             output << "cFetch(" << index << ")";
 4246                             produces = 0;
 4247                             break;
 4248                         }
 4249     #ifdef FP_SUPPORT_OPTIMIZER
 4250                         case cLog2by: n = "log2by"; params = 2; out_params = 1; break;
 4251                         case cPopNMov:
 4252                         {
 4253                             std::size_t a = ByteCode[++IP];
 4254                             std::size_t b = ByteCode[++IP];
 4255                             if(showExpression && b < stack.size())
 4256                             {
 4257                                 std::pair<int, std::string> stacktop(0, "?");
 4258                                 if(b < stack.size()) stacktop = stack[b];
 4259                                 stack.resize(a);
 4260                                 stack.push_back(stacktop);
 4261                             }
 4262                             output << "cPopNMov(" << a << ", " << b << ")";
 4263                             produces = 0;
 4264                             break;
 4265                         }
 4266                         case cNop:
 4267                             output << "nop"; params = 0; produces = 0;
 4268                             break;
 4269     #endif
 4270                         case cSinCos:
 4271                         {
 4272                             if(showExpression)
 4273                             {
 4274                                 std::pair<int, std::string> sin = stack.back();
 4275                                 std::pair<int, std::string> cos(
 4276                                     0, "cos(" + sin.second + ")");
 4277                                 sin.first = 0;
 4278                                 sin.second = "sin(" + sin.second + ")";
 4279                                 stack.back() = sin;
 4280                                 stack.push_back(cos);
 4281                             }
 4282                             output << "sincos";
 4283                             produces = 0;
 4284                             break;
 4285                         }
 4286                         case cSinhCosh:
 4287                         {
 4288                             if(showExpression)
 4289                             {
 4290                                 std::pair<int, std::string> sinh = stack.back();
 4291                                 std::pair<int, std::string> cosh(
 4292                                     0, "cosh(" + sinh.second + ")");
 4293                                 sinh.first = 0;
 4294                                 sinh.second = "sinh(" + sinh.second + ")";
 4295                                 stack.back() = sinh;
 4296                                 stack.push_back(cosh);
 4297                             }
 4298                             output << "sinhcosh";
 4299                             produces = 0;
 4300                             break;
 4301                         }
 4302                         case cAbsAnd: n = "abs_and"; break;
 4303                         case cAbsOr:  n = "abs_or"; break;
 4304                         case cAbsNot: n = "abs_not"; params = 1; break;
 4305                         case cAbsNotNot: n = "abs_notnot"; params = 1; break;
 4306                         case cDup:
 4307                         {
 4308                             if(showExpression)
 4309                                 stack.push_back(stack.back());
 4310                             output << "dup";
 4311                             produces = 0;
 4312                             break;
 4313                         }
 4314                         case cInv: n = "inv"; params = 1; break;
 4315                         case cSqr: n = "sqr"; params = 1; break;
 4316                         case cRDiv: n = "rdiv"; break;
 4317                         case cRSub: n = "rsub"; break;
 4318                         case cRSqrt: n = "rsqrt"; params = 1; break;
 4319 
 4320                         default:
 4321                             n = Functions[opcode-cAbs].name;
 4322                             params = Functions[opcode-cAbs].params;
 4323                             out_params = params != 1;
 4324                       }
 4325                   }
 4326             }
 4327         }
 4328         if(produces) output << n;
 4329         if(out_params) output << " (" << params << ")";
 4330         if(showExpression)
 4331         {
 4332             padLine(outputBuffer, 20);
 4333 
 4334             if(produces > 0)
 4335             {
 4336                 std::ostringstream buf;
 4337                 const char *paramsep = ",", *suff = "";
 4338                 int prio = 0; bool commutative = false;
 4339                 switch(opcode)
 4340                 {
 4341                   case cIf: buf << "if("; suff = ")";
 4342                       break;
 4343                   case cAbsIf: buf << "if("; suff = ")";
 4344                       break;
 4345                   case cOr:  prio = 6; paramsep = "|"; commutative = true;
 4346                       break;
 4347                   case cAnd: prio = 5; paramsep = "&"; commutative = true;
 4348                       break;
 4349                   case cAdd: prio = 4; paramsep = "+"; commutative = true;
 4350                       break;
 4351                   case cSub: prio = 4; paramsep = "-";
 4352                       break;
 4353                   case cMul: prio = 3; paramsep = "*"; commutative = true;
 4354                       break;
 4355                   case cDiv: prio = 3; paramsep = "/";
 4356                       break;
 4357                   case cPow: prio = 2; paramsep = "^";
 4358                       break;
 4359                   case cAbsOr:  prio = 6; paramsep = "|"; commutative = true;
 4360                       break;
 4361                   case cAbsAnd: prio = 5; paramsep = "&"; commutative = true;
 4362                       break;
 4363                   case cSqr: prio = 2; suff = "^2";
 4364                       break;
 4365                   case cNeg: buf << "(-("; suff = "))";
 4366                       break;
 4367                   case cNot: buf << "(!("; suff = "))";
 4368                       break;
 4369                   default: buf << n << '('; suff = ")";
 4370                 }
 4371 
 4372                 const char* sep = "";
 4373                 for(unsigned a=0; a<params; ++a)
 4374                 {
 4375                     buf << sep;
 4376                     if(stack.size() + a < params)
 4377                         buf << "?";
 4378                     else
 4379                     {
 4380                         const std::pair<int,std::string>& prev =
 4381                             stack[stack.size() - params + a];
 4382                         if(prio > 0 && (prev.first > prio ||
 4383                                         (prev.first==prio && !commutative)))
 4384                             buf << '(' << prev.second << ')';
 4385                         else
 4386                             buf << prev.second;
 4387                     }
 4388                     sep = paramsep;
 4389                 }
 4390                 if(stack.size() >= params)
 4391                     stack.resize(stack.size() - params);
 4392                 else
 4393                     stack.clear();
 4394                 buf << suff;
 4395                 stack.push_back(std::make_pair(prio, buf.str()));
 4396                 //if(n.size() <= 4 && !out_params) padLine(outputBuffer, 20);
 4397             }
 4398             //padLine(outputBuffer, 20);
 4399             output << "= ";
 4400             if(((opcode == cIf || opcode == cAbsIf) && params != 3)
 4401               || opcode == cJump
 4402     #ifdef FP_SUPPORT_OPTIMIZER
 4403               || opcode == cNop
 4404     #endif
 4405                 )
 4406                 output << "(void)";
 4407             else if(stack.empty())
 4408                 output << "[?] ?";
 4409             else
 4410                 output << '[' << (stack.size()-1) << ']'
 4411                        << stack.back().second;
 4412         }
 4413 
 4414         if(showExpression)
 4415         {
 4416             dest << outputBuffer.str() << std::endl;
 4417             outputBuffer.str("");
 4418         }
 4419         else
 4420             output << std::endl;
 4421     }
 4422     dest << std::flush;
 4423 }
 4424 #endif
 4425 
 4426 
 4427 #ifndef FP_SUPPORT_OPTIMIZER
 4428 template<typename Value_t>
 4429 void FunctionParserBase<Value_t>::Optimize()
 4430 {
 4431     // Do nothing if no optimizations are supported.
 4432 }
 4433 #endif
 4434 
 4435 
 4436 #define FUNCTIONPARSER_INSTANTIATE_CLASS(type) \
 4437     template class FunctionParserBase< type >;
 4438 
 4439 #ifndef FP_DISABLE_DOUBLE_TYPE
 4440 FUNCTIONPARSER_INSTANTIATE_CLASS(double)
 4441 #endif
 4442 
 4443 #ifdef FP_SUPPORT_FLOAT_TYPE
 4444 FUNCTIONPARSER_INSTANTIATE_CLASS(float)
 4445 #endif
 4446 
 4447 #ifdef FP_SUPPORT_LONG_DOUBLE_TYPE
 4448 FUNCTIONPARSER_INSTANTIATE_CLASS(long double)
 4449 #endif
 4450 
 4451 #ifdef FP_SUPPORT_LONG_INT_TYPE
 4452 FUNCTIONPARSER_INSTANTIATE_CLASS(long)
 4453 #endif
 4454 
 4455 #ifdef FP_SUPPORT_MPFR_FLOAT_TYPE
 4456 FUNCTIONPARSER_INSTANTIATE_CLASS(MpfrFloat)
 4457 #endif
 4458 
 4459 #ifdef FP_SUPPORT_GMP_INT_TYPE
 4460 FUNCTIONPARSER_INSTANTIATE_CLASS(GmpInt)
 4461 #endif
 4462 
 4463 #ifdef FP_SUPPORT_COMPLEX_DOUBLE_TYPE
 4464 FUNCTIONPARSER_INSTANTIATE_CLASS(std::complex<double>)
 4465 #endif
 4466 
 4467 #ifdef FP_SUPPORT_COMPLEX_FLOAT_TYPE
 4468 FUNCTIONPARSER_INSTANTIATE_CLASS(std::complex<float>)
 4469 #endif
 4470 
 4471 #ifdef FP_SUPPORT_COMPLEX_LONG_DOUBLE_TYPE
 4472 FUNCTIONPARSER_INSTANTIATE_CLASS(std::complex<long double>)
 4473 #endif