diff --git a/ADOL-C/boost-test/traceOperatorScalar.cpp b/ADOL-C/boost-test/traceOperatorScalar.cpp index f58e4ff5..0136bfda 100644 --- a/ADOL-C/boost-test/traceOperatorScalar.cpp +++ b/ADOL-C/boost-test/traceOperatorScalar.cpp @@ -4514,6 +4514,99 @@ BOOST_AUTO_TEST_CASE(ErfOperator_FOS_Reverse) myfree1(z); } +BOOST_AUTO_TEST_CASE(ErfcOperator_ZOS_Forward) +{ + double a = 7.1, aout; + adouble ad; + + trace_on(1); + ad <<= a; + + ad = erfc(ad); + + ad >>= aout; + trace_off(); + + a = std::erfc(a); + + BOOST_TEST(aout == a, tt::tolerance(tol)); + + double *x = myalloc1(1); + double *y = myalloc1(1); + + *x = 7.1; + + zos_forward(1, 1, 1, 0, x, y); + + BOOST_TEST(*y == a, tt::tolerance(tol)); + + myfree1(x); + myfree1(y); +} + +BOOST_AUTO_TEST_CASE(ErfcOperator_FOS_Forward) +{ + double a = 7.1, aout; + adouble ad; + + trace_on(1); + ad <<= a; + + ad = erfc(ad); + + ad >>= aout; + trace_off(); + + double aDerivative = -2. / std::sqrt(std::acos(-1.)) * std::exp(-a * a); + a = std::erfc(a); + + double *x = myalloc1(1); + double *xd = myalloc1(1); + double *y = myalloc1(1); + double *yd = myalloc1(1); + + *x = 7.1; + *xd = 1.; + + fos_forward(1, 1, 1, 0, x, xd, y, yd); + + BOOST_TEST(*y == a, tt::tolerance(tol)); + BOOST_TEST(*yd == aDerivative, tt::tolerance(tol)); + + myfree1(x); + myfree1(xd); + myfree1(y); + myfree1(yd); +} + +BOOST_AUTO_TEST_CASE(ErfcOperator_FOS_Reverse) +{ + double a = 7.1, aout; + adouble ad; + + trace_on(1, 1); + ad <<= a; + + ad = erfc(ad); + + ad >>= aout; + trace_off(); + + double aDerivative = -2. / std::sqrt(std::acos(-1.)) * std::exp(-a * a); + + double *u = myalloc1(1); + double *z = myalloc1(1); + + *u = 1.; + + fos_reverse(1, 1, 1, u, z); + + BOOST_TEST(*z == aDerivative, tt::tolerance(tol)); + + myfree1(u); + myfree1(z); +} + BOOST_AUTO_TEST_CASE(EqOperator_ZOS_Forward) { double a = 10.01, aout; diff --git a/ADOL-C/boost-test/traceOperatorVector.cpp b/ADOL-C/boost-test/traceOperatorVector.cpp index 3867ba0a..0fb6c594 100644 --- a/ADOL-C/boost-test/traceOperatorVector.cpp +++ b/ADOL-C/boost-test/traceOperatorVector.cpp @@ -3394,6 +3394,75 @@ BOOST_AUTO_TEST_CASE(ErfOperator_FOV_Reverse) myfree2(z); } +BOOST_AUTO_TEST_CASE(ErfcOperator_FOV_Forward) +{ + double a = 7.1, aout; + adouble ad; + + trace_on(1); + ad <<= a; + + ad = erfc(ad); + + ad >>= aout; + trace_off(); + + double aDerivative = -2. / std::sqrt(std::acos(-1.)) * std::exp(-a*a); + a = std::erfc(a); + + double *x = myalloc1(1); + double **xd = myalloc2(1, 2); + double *y = myalloc1(1); + double **yd = myalloc2(1, 2); + + x[0] = 7.1; + + for (int i = 0; i < 2; i++) { + xd[0][i] = std::pow(3., i*2.); + } + + fov_forward(1, 1, 1, 2, x, xd, y, yd); + + BOOST_TEST(*y == a, tt::tolerance(tol)); + BOOST_TEST(yd[0][0] == aDerivative, tt::tolerance(tol)); + BOOST_TEST(yd[0][1] == aDerivative * std::pow(3., 2.), tt::tolerance(tol)); + + myfree1(x); + myfree2(xd); + myfree1(y); + myfree2(yd); +} + +BOOST_AUTO_TEST_CASE(ErfcOperator_FOV_Reverse) +{ + double a = 7.1, aout; + adouble ad; + + trace_on(1, 1); + ad <<= a; + + ad = erfc(ad); + + ad >>= aout; + trace_off(); + + double aDerivative = -2. / std::sqrt(std::acos(-1.)) * std::exp(-a*a); + + double **u = myalloc2(2, 1); + double **z = myalloc2(2, 1); + + u[0][0] = 1.; + u[1][0] = -1.1; + + fov_reverse(1, 1, 1, 2, u, z); + + BOOST_TEST(z[0][0] == aDerivative, tt::tolerance(tol)); + BOOST_TEST(z[1][0] == -1.1*aDerivative, tt::tolerance(tol)); + + myfree2(u); + myfree2(z); +} + BOOST_AUTO_TEST_CASE(EqPlusOperator_FOV_Forward) { double a = 5.132, aout; diff --git a/ADOL-C/boost-test/tracelessCompositeTests.cpp b/ADOL-C/boost-test/tracelessCompositeTests.cpp index 1ad29cdc..47172f70 100644 --- a/ADOL-C/boost-test/tracelessCompositeTests.cpp +++ b/ADOL-C/boost-test/tracelessCompositeTests.cpp @@ -282,6 +282,32 @@ BOOST_AUTO_TEST_CASE(CompositeErfFabs_Traceless) BOOST_TEST(ax1.getADValue(1) == x2Derivative, tt::tolerance(tol)); } +BOOST_AUTO_TEST_CASE(CompositeErfcFabs_Traceless) +{ + double x1 = 4.56, x2 = 6.45; + adouble ax1 = x1, ax2 = x2; + + ax1.setADValue(0, 1.); + ax2.setADValue(1, 1.); + + ax1 = adtl::erfc(adtl::fabs(ax1 - 2.)*adtl::sinh(ax2)); + + double x1Derivative = -2. / std::sqrt(std::acos(-1.)) + * std::exp(-std::pow(std::fabs(x1 - 2.) + * std::sinh(x2), 2)) + * std::sinh(x2); + double x2Derivative = -2. / std::sqrt(std::acos(-1.)) + * std::exp(-std::pow(std::fabs(x1 - 2.) + * std::sinh(x2), 2)) + * std::fabs(x1 - 2.) * std::cosh(x2); + + x1 = std::erfc(std::fabs(x1 - 2.)*std::sinh(x2)); + + BOOST_TEST(ax1.getValue() == x1, tt::tolerance(tol)); + BOOST_TEST(ax1.getADValue(0) == x1Derivative, tt::tolerance(tol)); + BOOST_TEST(ax1.getADValue(1) == x2Derivative, tt::tolerance(tol)); +} + BOOST_AUTO_TEST_CASE(ExpTrigSqrt_Traceless) { double x1 = 1.2, x2 = 2.1; diff --git a/ADOL-C/boost-test/tracelessOperatorScalar.cpp b/ADOL-C/boost-test/tracelessOperatorScalar.cpp index b543f6e3..00737902 100644 --- a/ADOL-C/boost-test/tracelessOperatorScalar.cpp +++ b/ADOL-C/boost-test/tracelessOperatorScalar.cpp @@ -1436,6 +1436,36 @@ BOOST_AUTO_TEST_CASE(ErfOperatorDerivative) BOOST_TEST(ad.getADValue(0) == aDerivative, tt::tolerance(tol)); } +/* The complementary error function erfc(a) is defined as + * 1.0 - erf(a). + */ +BOOST_AUTO_TEST_CASE(ErfcOperatorPrimal) +{ + double a = 7.1; + adouble ad = a; + + a = std::erfc(a); + ad = adtl::erfc(ad); + + BOOST_TEST(ad.getValue() == a, tt::tolerance(tol)); +} + +BOOST_AUTO_TEST_CASE(ErfcOperatorDerivative) +{ + double a = 7.1; + adouble ad = a; + + /* The inverse cosine is used for pi so that no additional math library + * import is necessary, see also the implementation in adtl.h. + */ + double aDerivative = -2. / std::sqrt(std::acos(-1.)) * std::exp(-a * a); + + ad.setADValue(0, 1.); + ad = adtl::erfc(ad); + + BOOST_TEST(ad.getADValue(0) == aDerivative, tt::tolerance(tol)); +} + /* Test the primitive non-temporary operations =, +=, -=, *=, /=. */ BOOST_AUTO_TEST_CASE(EqOperatorPrimal_Derivative) diff --git a/ADOL-C/boost-test/tracelessOperatorVector.cpp b/ADOL-C/boost-test/tracelessOperatorVector.cpp index 24d4d7fd..e060098b 100644 --- a/ADOL-C/boost-test/tracelessOperatorVector.cpp +++ b/ADOL-C/boost-test/tracelessOperatorVector.cpp @@ -868,6 +868,21 @@ BOOST_AUTO_TEST_CASE(ErfOperatorDerivativeVectorMode) BOOST_TEST(ad.getADValue(1) == aDerivative*(-2.5), tt::tolerance(tol)); } +BOOST_AUTO_TEST_CASE(ErfcOperatorDerivativeVectorMode) +{ + double a = 7.1; + adouble ad = a; + + double aDerivative = -2. / std::sqrt(std::acos(-1.)) * std::exp(-a * a); + + ad.setADValue(0, 1.); + ad.setADValue(1, -2.5); + ad = adtl::erfc(ad); + + BOOST_TEST(ad.getADValue(0) == aDerivative, tt::tolerance(tol)); + BOOST_TEST(ad.getADValue(1) == aDerivative*(-2.5), tt::tolerance(tol)); +} + BOOST_AUTO_TEST_CASE(EqOperatorDerivativeVectorMode) { double a = 10.01; diff --git a/ADOL-C/c_interface/ADOLC_TB_interface.cpp b/ADOL-C/c_interface/ADOLC_TB_interface.cpp index 1e7036d4..6dfcea24 100644 --- a/ADOL-C/c_interface/ADOLC_TB_interface.cpp +++ b/ADOL-C/c_interface/ADOLC_TB_interface.cpp @@ -261,6 +261,10 @@ extern "C" { return new adouble(erf(*static_cast(a))); } + TBAdoubleHandle tb_erfc(TBAdoubleHandle a) + { + return new adouble(erfc(*static_cast(a))); + } } /* diff --git a/ADOL-C/c_interface/ADOLC_TB_interface.h b/ADOL-C/c_interface/ADOLC_TB_interface.h index 438d42e4..bce1d0ef 100644 --- a/ADOL-C/c_interface/ADOLC_TB_interface.h +++ b/ADOL-C/c_interface/ADOLC_TB_interface.h @@ -115,6 +115,7 @@ extern "C" TBAdoubleHandle tb_floor(TBAdoubleHandle a); TBAdoubleHandle tb_ldexp(TBAdoubleHandle a, int n); TBAdoubleHandle tb_erf(TBAdoubleHandle a); + TBAdoubleHandle tb_erfc(TBAdoubleHandle a); #ifdef __cplusplus } #endif @@ -139,4 +140,4 @@ extern "C" } #endif -#endif // ADOLC_TB_INTERFACE_H \ No newline at end of file +#endif // ADOLC_TB_INTERFACE_H diff --git a/ADOL-C/c_interface/ADOLC_TL_interface.cpp b/ADOL-C/c_interface/ADOLC_TL_interface.cpp index de85c659..05beab39 100644 --- a/ADOL-C/c_interface/ADOLC_TL_interface.cpp +++ b/ADOL-C/c_interface/ADOLC_TL_interface.cpp @@ -292,4 +292,8 @@ extern "C" { return new adtl::adouble(erf(*static_cast(a))); } + TLAdoubleHandle tl_erfc(TLAdoubleHandle a) + { + return new adtl::adouble(erfc(*static_cast(a))); + } } diff --git a/ADOL-C/c_interface/ADOLC_TL_interface.h b/ADOL-C/c_interface/ADOLC_TL_interface.h index 5558463e..e329644b 100644 --- a/ADOL-C/c_interface/ADOLC_TL_interface.h +++ b/ADOL-C/c_interface/ADOLC_TL_interface.h @@ -122,8 +122,9 @@ extern "C" TLAdoubleHandle tl_floor(TLAdoubleHandle a); TLAdoubleHandle tl_ldexp(TLAdoubleHandle a, const int n); TLAdoubleHandle tl_erf(TLAdoubleHandle a); + TLAdoubleHandle tl_erfc(TLAdoubleHandle a); #ifdef __cplusplus } #endif -#endif // ADOLC_TL_INTERFACE_H \ No newline at end of file +#endif // ADOLC_TL_INTERFACE_H diff --git a/ADOL-C/include/adolc/adoublecuda.h b/ADOL-C/include/adolc/adoublecuda.h index 78dc15b8..32acc231 100644 --- a/ADOL-C/include/adolc/adoublecuda.h +++ b/ADOL-C/include/adolc/adoublecuda.h @@ -152,6 +152,7 @@ class adouble { CUDADEV friend adouble ldexp (const double v, const adouble &a); CUDADEV friend double frexp (const adouble &a, int* v); CUDADEV friend adouble erf (const adouble &a); + CUDADEV friend adouble erfc (const adouble &a); /******************* nontemporary results ***************************/ @@ -794,6 +795,16 @@ CUDADEV adouble erf (const adouble &a) { return tmp; } +CUDADEV adouble erfc (const adouble &a) { + adouble tmp; + tmp.val=ADOLC_MATH_NSP_ERF::erfc(a.val); + double tmp2 = -2.0 / + ADOLC_MATH_NSP_ERF::sqrt(ADOLC_MATH_NSP::acos(-1.0)) * + ADOLC_MATH_NSP_ERF::exp(-a.val*a.val); + FOR_I_EQ_0_LT_NUMDIR + tmp.ADVAL_I=tmp2*a.ADVAL_I; + return tmp; +} /******************* nontemporary results *********************************/ CUDADEV void adouble::operator = (const double v) { diff --git a/ADOL-C/include/adolc/adtl.h b/ADOL-C/include/adolc/adtl.h index c62cde47..88865abb 100644 --- a/ADOL-C/include/adolc/adtl.h +++ b/ADOL-C/include/adolc/adtl.h @@ -137,6 +137,7 @@ class adouble { inline friend adouble ldexp (const double v, const adouble &a); inline friend double frexp (const adouble &a, int* v); inline friend adouble erf (const adouble &a); + inline friend adouble erfc (const adouble &a); inline friend void condassign( adouble &res, const adouble &cond, const adouble &arg1, const adouble &arg2 ); @@ -866,6 +867,17 @@ inline adouble erf (const adouble &a) { return tmp; } +inline adouble erfc (const adouble &a) { + adouble tmp; + tmp.PRIMAL_VALUE=ADOLC_MATH_NSP_ERF::erfc(a.PRIMAL_VALUE); + double tmp2 = -2.0 / + ADOLC_MATH_NSP_ERF::sqrt(ADOLC_MATH_NSP::acos(-1.0)) * + ADOLC_MATH_NSP_ERF::exp(-a.PRIMAL_VALUE*a.PRIMAL_VALUE); + FOR_I_EQ_1_LTEQ_NUMDIR + tmp.ADVAL_I=tmp2*a.ADVAL_I; + return tmp; +} + inline void condassign( adouble &res, const adouble &cond, const adouble &arg1, const adouble &arg2 ) { if (cond.getValue() > 0) diff --git a/ADOL-C/include/adolc/adtl_hov.h b/ADOL-C/include/adolc/adtl_hov.h index 44533df3..b416bdb2 100644 --- a/ADOL-C/include/adolc/adtl_hov.h +++ b/ADOL-C/include/adolc/adtl_hov.h @@ -147,6 +147,7 @@ class adouble { inline friend adouble ldexp (const double v, const adouble &a); inline friend double frexp (const adouble &a, int* v); inline friend adouble erf (const adouble &a); + inline friend adouble erfc (const adouble &a); inline friend void condassign( adouble &res, const adouble &cond, const adouble &arg1, const adouble &arg2 ); @@ -1937,6 +1938,26 @@ inline adouble erf (const adouble &a) { return tmp; } +inline adouble erfc (const adouble &a) { + adouble tmp; + if (unlikely(!adouble::_do_val() && adouble::_do_adval())) { + fprintf(DIAG_OUT, "ADOL-C error: Traceless: Incorrect mode, call setMode(enum Mode mode)\n"); + throw logic_error("incorrect function call, errorcode=1"); + } + if (do_val()) + tmp.val=ADOLC_MATH_NSP_ERF::erfc(a.val); + if (likely(adouble::_do_adval() && adouble::_do_val())) { + double tmp2 = -2.0 / + ADOLC_MATH_NSP_ERF::sqrt(ADOLC_MATH_NSP::acos(-1.0)) * + ADOLC_MATH_NSP_ERF::exp(-a.val*a.val); + FOR_I_EQ_0_LT_NUMDIR + tmp.ADVAL_I=tmp2*a.ADVAL_I; + } + if (do_indo()) + tmp.add_to_pattern( a.get_pattern() ) ; + return tmp; +} + inline void condassign( adouble &res, const adouble &cond, const adouble &arg1, const adouble &arg2 ) { if (no_do_val()) { diff --git a/ADOL-C/include/adolc/adtl_indo.h b/ADOL-C/include/adolc/adtl_indo.h index 349ec7cd..f9444268 100644 --- a/ADOL-C/include/adolc/adtl_indo.h +++ b/ADOL-C/include/adolc/adtl_indo.h @@ -129,6 +129,7 @@ class adouble { inline friend adouble ldexp (const double v, const adouble &a); inline friend double frexp (const adouble &a, int* v); inline friend adouble erf (const adouble &a); + inline friend adouble erfc (const adouble &a); inline friend void condassign( adouble &res, const adouble &cond, const adouble &arg1, const adouble &arg2 ); @@ -690,6 +691,13 @@ inline adouble erf (const adouble &a) { return tmp; } +inline adouble erfc (const adouble &a) { + adouble tmp; + tmp.val=ADOLC_MATH_NSP_ERF::erfc(a.val); + tmp.add_to_pattern( a.get_pattern() ) ; + return tmp; +} + inline void condassign( adouble &res, const adouble &cond, const adouble &arg1, const adouble &arg2 ) { if (cond.getValue() > 0) diff --git a/ADOL-C/include/adolc/internal/adubfunc.h b/ADOL-C/include/adolc/internal/adubfunc.h index 02315d57..6650892e 100644 --- a/ADOL-C/include/adolc/internal/adubfunc.h +++ b/ADOL-C/include/adolc/internal/adubfunc.h @@ -72,6 +72,7 @@ friend ADOLC_DLL_EXPORT adub acosh ( const badouble& ); friend ADOLC_DLL_EXPORT adub atanh ( const badouble& ); friend ADOLC_DLL_EXPORT adub erf ( const badouble& ); + friend ADOLC_DLL_EXPORT adub erfc ( const badouble& ); friend ADOLC_DLL_EXPORT adub fabs ( const badouble& ); friend ADOLC_DLL_EXPORT adub ceil ( const badouble& ); diff --git a/ADOL-C/include/adolc/internal/paramfunc.h b/ADOL-C/include/adolc/internal/paramfunc.h index d1373e4e..40cc84f2 100644 --- a/ADOL-C/include/adolc/internal/paramfunc.h +++ b/ADOL-C/include/adolc/internal/paramfunc.h @@ -96,6 +96,7 @@ inline friend ADOLC_DLL_EXPORT adub acosh ( const pdouble& ); inline friend ADOLC_DLL_EXPORT adub atanh ( const pdouble& ); inline friend ADOLC_DLL_EXPORT adub erf ( const pdouble& ); + inline friend ADOLC_DLL_EXPORT adub erfc ( const pdouble& ); inline friend ADOLC_DLL_EXPORT adub fabs ( const pdouble& ); inline friend ADOLC_DLL_EXPORT adub ceil ( const pdouble& ); diff --git a/ADOL-C/include/adolc/param.h b/ADOL-C/include/adolc/param.h index 4ff1faed..27a19768 100644 --- a/ADOL-C/include/adolc/param.h +++ b/ADOL-C/include/adolc/param.h @@ -185,6 +185,7 @@ inline adub asinh ( const pdouble& p) { return asinh(adub(p)); } inline adub acosh ( const pdouble& p) { return acosh(adub(p)); } inline adub atanh ( const pdouble& p) { return atanh(adub(p)); } inline adub erf ( const pdouble& p) { return erf(adub(p)); } +inline adub erfc ( const pdouble& p) { return erfc(adub(p)); } inline adub fabs ( const pdouble& p) { return fabs(adub(p)); } inline adub ceil ( const pdouble& p) { return ceil(adub(p)); } inline adub floor ( const pdouble& p) { return floor(adub(p)); } diff --git a/ADOL-C/src/adouble.cpp b/ADOL-C/src/adouble.cpp index 7107edc3..dda45571 100644 --- a/ADOL-C/src/adouble.cpp +++ b/ADOL-C/src/adouble.cpp @@ -2999,6 +2999,56 @@ adub erf( const badouble& x ) { return locat; } +adub erfc( const badouble& x ) { + ADOLC_OPENMP_THREAD_NUMBER; + ADOLC_OPENMP_GET_THREAD_NUMBER; + locint locat = next_loc(); + double coval = ADOLC_MATH_NSP_ERF::erfc(ADOLC_GLOBAL_TAPE_VARS.store[x.loc()]); + + adouble y = -2.0 / + ADOLC_MATH_NSP_ERF::sqrt(ADOLC_MATH_NSP::acos(-1.0))*exp(-x*x); + + if (ADOLC_CURRENT_TAPE_INFOS.traceFlag) { // old: write_quad(erfc_op,locat,x.loc(),y.loc()); +#if defined(ADOLC_TRACK_ACTIVITY) + if (ADOLC_GLOBAL_TAPE_VARS.actStore[x.loc()]) { // y will have same activity as x and can be considered as second input here +#endif + put_op(erfc_op); + ADOLC_PUT_LOCINT(x.loc()); // = arg1 + ADOLC_PUT_LOCINT(y.loc()); // = arg2 + ADOLC_PUT_LOCINT(locat); // = res + + ++ADOLC_CURRENT_TAPE_INFOS.numTays_Tape; + if (ADOLC_CURRENT_TAPE_INFOS.keepTaylors) + ADOLC_WRITE_SCAYLOR(ADOLC_GLOBAL_TAPE_VARS.store[locat]); +#if defined(ADOLC_TRACK_ACTIVITY) + } else if (ADOLC_GLOBAL_TAPE_VARS.actStore[locat]) { + if (coval == 0.0) { + put_op(assign_d_zero); + ADOLC_PUT_LOCINT(locat); + } else if (coval == 1.0) { + put_op(assign_d_one); + ADOLC_PUT_LOCINT(locat); + } else { + put_op(assign_d); + ADOLC_PUT_LOCINT(locat); + ADOLC_PUT_VAL(coval); + } + + ++ADOLC_CURRENT_TAPE_INFOS.numTays_Tape; + if (ADOLC_CURRENT_TAPE_INFOS.keepTaylors) + ADOLC_WRITE_SCAYLOR(ADOLC_GLOBAL_TAPE_VARS.store[locat]); + } +#endif + } + + ADOLC_GLOBAL_TAPE_VARS.store[locat] = coval; +#if defined(ADOLC_TRACK_ACTIVITY) + ADOLC_GLOBAL_TAPE_VARS.actStore[locat] = ADOLC_GLOBAL_TAPE_VARS.actStore[x.loc()]; +#endif + ADOLC_OPENMP_RESTORE_THREAD_NUMBER; + return locat; +} + /*--------------------------------------------------------------------------*/ /* Fabs Function (NOTE: This function is also nondifferentiable at x=0) */ adub fabs ( const badouble& x ) { diff --git a/ADOL-C/src/fo_rev.c b/ADOL-C/src/fo_rev.c index db8c1dea..6fd04b06 100644 --- a/ADOL-C/src/fo_rev.c +++ b/ADOL-C/src/fo_rev.c @@ -1484,6 +1484,7 @@ int int_reverse_safe( case acosh_op: /* acosh_op */ case atanh_op: /* atanh_op */ case erf_op: /* erf_op */ + case erfc_op: /* erfc_op */ res = get_locint_r(); arg2 = get_locint_r(); arg1 = get_locint_r(); diff --git a/ADOL-C/src/ho_rev.c b/ADOL-C/src/ho_rev.c index 779e1d79..de1404e3 100644 --- a/ADOL-C/src/ho_rev.c +++ b/ADOL-C/src/ho_rev.c @@ -1635,6 +1635,7 @@ int hov_ti_reverse( case acosh_op: /* acosh_op */ case atanh_op: /* atanh_op */ case erf_op: /* erf_op */ + case erfc_op: /* erfc_op */ res = get_locint_r(); arg2 = get_locint_r(); arg1 = get_locint_r(); diff --git a/ADOL-C/src/oplate.h b/ADOL-C/src/oplate.h index 478225fe..6d518629 100644 --- a/ADOL-C/src/oplate.h +++ b/ADOL-C/src/oplate.h @@ -21,7 +21,7 @@ /* opcodes */ enum OPCODES { - death_not = 0, + death_not, assign_ind, assign_dep, assign_a, @@ -40,7 +40,7 @@ enum OPCODES { mult_d_a, div_a_a, div_d_a, - exp_op = 19, + exp_op, cos_op, sin_op, atan_op, @@ -52,13 +52,13 @@ enum OPCODES { asinh_op, acosh_op, atanh_op, - gen_quad = 31, + gen_quad, end_of_tape, start_of_tape, end_of_op, end_of_int, end_of_val, - cond_assign = 37, + cond_assign, cond_assign_s, take_stock_op, assign_d_one, @@ -67,7 +67,7 @@ enum OPCODES { decr_a, neg_sign_a, pos_sign_a, - min_op = 46, + min_op, abs_val, eq_zero, neq_zero, @@ -75,18 +75,19 @@ enum OPCODES { gt_zero, ge_zero, lt_zero, - eq_plus_prod = 54, + eq_plus_prod, eq_min_prod, erf_op, + erfc_op, ceil_op, floor_op, - ext_diff = 59, + ext_diff, ext_diff_iArr, ignore_me, ext_diff_v2, - cond_eq_assign = 63, + cond_eq_assign, cond_eq_assign_s, - subscript = 80, + subscript, subscript_ref, ref_assign_d_zero, ref_assign_d_one, @@ -104,7 +105,7 @@ enum OPCODES { ref_copyout, ref_cond_assign, ref_cond_assign_s, - assign_p = 98, + assign_p, eq_plus_p, eq_min_p, eq_mult_p, @@ -122,21 +123,21 @@ enum OPCODES { vec_copy, vec_dot, vec_axpy, - ref_cond_eq_assign = 116, + ref_cond_eq_assign, ref_cond_eq_assign_s, - eq_a_p = 119, + eq_a_p, neq_a_p, le_a_p, gt_a_p, ge_a_p, lt_a_p, - eq_a_a = 125, + eq_a_a, neq_a_a, le_a_a, gt_a_a, ge_a_a, lt_a_a, - ampi_send = 131, + ampi_send, ampi_recv, ampi_isend, ampi_irecv, @@ -152,7 +153,7 @@ enum OPCODES { ampi_reduce, ampi_allreduce, medi_call, - cbrt_op + cbrt_op, }; /****************************************************************************/ diff --git a/ADOL-C/src/tapedoc/tapedoc.c b/ADOL-C/src/tapedoc/tapedoc.c index 26300c82..c2abd59f 100644 --- a/ADOL-C/src/tapedoc/tapedoc.c +++ b/ADOL-C/src/tapedoc/tapedoc.c @@ -968,6 +968,25 @@ void tape_doc(short tnum, /* tape id */ filewrite(operation,"erf op",3,loc_a,val_a,0,cst_d); break; + /*--------------------------------------------------------------------------*/ + case erfc_op: /* erfc_op */ + arg1 = get_locint_f(); + arg2 = get_locint_f(); + res = get_locint_f(); + loc_a[0]=arg1; + loc_a[1]=arg2; + loc_a[2]=res; +#ifdef ADOLC_TAPE_DOC_VALUES + val_a[0]=dp_T0[arg1]; + dp_T0[res] = erfc(dp_T0[arg1]); + ADOLC_OPENMP_RESTORE_THREAD_NUMBER; + val_a[1]=dp_T0[arg2]; + val_a[2]=dp_T0[res]; +#endif + filewrite(operation,"erfc op",3,loc_a,val_a,0,cst_d); + break; + + /*--------------------------------------------------------------------------*/ case log_op: /* log_op */ arg = get_locint_f(); diff --git a/ADOL-C/src/uni5_for.c b/ADOL-C/src/uni5_for.c index 34ef864e..9becf479 100644 --- a/ADOL-C/src/uni5_for.c +++ b/ADOL-C/src/uni5_for.c @@ -3521,6 +3521,66 @@ int hov_forward( TRES_FOINC = dp_T0[arg2] * TARG1_INC; +#if defined(_HIGHER_ORDER_) + Targ2OP = Targ2; + + *Tres *= (i+1); + for (int j=0;j