From 6da3fc5a89b60d1f01f885f4c558e7f42ed7fe30 Mon Sep 17 00:00:00 2001 From: Neil Lawrence Date: Wed, 27 Nov 2013 11:17:33 +0000 Subject: [PATCH] Added gradient of sympy kernel, seems to pass tests, but know it's not numerically stable. Checking in before making numerically stable. --- GPy/kern/parts/sympy_helpers.cpp | 20 ++++++++++++++++++-- 1 file changed, 18 insertions(+), 2 deletions(-) diff --git a/GPy/kern/parts/sympy_helpers.cpp b/GPy/kern/parts/sympy_helpers.cpp index 9f30eea9..9b0d5885 100644 --- a/GPy/kern/parts/sympy_helpers.cpp +++ b/GPy/kern/parts/sympy_helpers.cpp @@ -170,9 +170,25 @@ double dh_dl(double t, double tprime, double d_i, double d_j, double l){ } double dh_dt(double t, double tprime, double d_i, double d_j, double l){ - return 0.0; + // compute gradient of h function with respect to t. + double diff_t = t - tprime; + double half_l_di = 0.5*l*d_i; + double arg_1 = half_l_di + tprime/l; + double arg_2 = half_l_di - diff_t/l; + double ln_part_1 = ln_diff_erf(arg_1, arg_2); + arg_2 = half_l_di - t/l; + double ln_part_2 = ln_diff_erf(half_l_di, arg_2); + + return (d_i*(erf(d_i*l/2) - erf(d_i*l/2 - t/l))*exp(-d_i*t - d_j*tprime) - d_i*(erf(d_i*l/2 + tprime/l) - erf(d_i*l/2 - (t - tprime)/l))*exp(-d_i*(t - tprime)) + 2*exp(-d_i*(t - tprime) - pow(d_i*l/2 - (t - tprime)/l, 2))/(sqrt(M_PI)*l) - 2*exp(-d_i*t - d_j*tprime - pow(d_i*l/2 - t/l,2))/(sqrt(M_PI)*l))*exp(d_i*l/2*d_i*l/2)/(d_i + d_j); } double dh_dtprime(double t, double tprime, double d_i, double d_j, double l){ - return 0.0; + // compute gradient of h function with respect to tprime. + double diff_t = t - tprime; + double half_l_di = 0.5*l*d_i; + double arg_1 = half_l_di + tprime/l; + double arg_2 = half_l_di - diff_t/l; + double ln_part_1 = ln_diff_erf(arg_1, arg_2); + + return (d_i*(erf(d_i*l/2 + tprime/l) - erf(d_i*l/2 - (t - tprime)/l))*exp(-d_i*(t - tprime)) + d_j*(erf(d_i*l/2) - erf(d_i*l/2 - t/l))*exp(-d_i*t - d_j*tprime) + (-2*exp(-pow(d_i*l/2 - (t - tprime)/l,2)) + 2*exp(-pow(d_i*l/2 + tprime/l,2)))*exp(-d_i*(t - tprime))/(sqrt(M_PI)*l))*exp(d_i*l/2*d_i*l/2)/(d_i + d_j); }