diff --git a/applications/fastsum/fastsum.m b/applications/fastsum/fastsum.m
index 2358b985..3aefc435 100644
--- a/applications/fastsum/fastsum.m
+++ b/applications/fastsum/fastsum.m
@@ -48,6 +48,7 @@
 %   'cot'                     K(x) = cot(cx)
 %   'one_over_cube'           K(x) = 1/x^3
 %   'laplacian_rbf'           K(x) = EXP(-|x|/c)
+%   'der_laplacian_rbf'       K(x) = |x|/c EXP(-|x|/c)
 %   'xx_gaussian'             K(x) = x^2/c^2 EXP(-x^2/c^2)
 %   'absx'                    K(x) = |x|
 %
diff --git a/applications/fastsum/fastsum_matlab.c b/applications/fastsum/fastsum_matlab.c
index 365eab23..952cb91f 100644
--- a/applications/fastsum/fastsum_matlab.c
+++ b/applications/fastsum/fastsum_matlab.c
@@ -122,6 +122,8 @@ int main(int argc, char **argv)
       kernel = log_sin;
     else if (strcmp(s, "laplacian_rbf") == 0)
       kernel = laplacian_rbf;
+    else if (strcmp(s, "der_laplacian_rbf") == 0)
+      kernel = der_laplacian_rbf;
     else if (strcmp(s, "xx_gaussian") == 0)
       kernel = xx_gaussian;
     else if (strcmp(s, "absx") == 0)
diff --git a/applications/fastsum/fastsum_test.c b/applications/fastsum/fastsum_test.c
index 016e6eb4..a9a0ea99 100644
--- a/applications/fastsum/fastsum_test.c
+++ b/applications/fastsum/fastsum_test.c
@@ -124,6 +124,8 @@ int main(int argc, char **argv)
       kernel = log_sin;
     else if (strcmp(s, "laplacian_rbf") == 0)
       kernel = laplacian_rbf;
+    else if (strcmp(s, "der_laplacian_rbf") == 0)
+      kernel = der_laplacian_rbf;
     else if (strcmp(s, "xx_gaussian") == 0)
       kernel = xx_gaussian;
     else if (strcmp(s, "absx") == 0)
diff --git a/applications/fastsum/kernels.c b/applications/fastsum/kernels.c
index f17f66b3..de27b90a 100644
--- a/applications/fastsum/kernels.c
+++ b/applications/fastsum/kernels.c
@@ -431,6 +431,22 @@ C laplacian_rbf(R x, int der, const R *param)    /* K(x)=EXP(-|x|/c) */
   return value;
 }
 
+C der_laplacian_rbf(R x, int der, const R *param)    /* K(x)=|x|/c EXP(-|x|/c) */
+{
+  R c = param[0];
+  R value = K(0.0);
+
+  switch (der)
+  {
+    case  0 : value = (FABS(x)/c)*EXP(-FABS(x)/c); break;
+    default:
+        value = (POW(K(-1.0),(R)der))*((FABS(x)-(R)der*c)/POW(c,(R)der+1))*EXP(-FABS(x)/c);
+        value *= 1 - 2 * ((x < K(0.0)) && (der % 2));
+  }
+
+  return value;
+}
+
 C xx_gaussian(R x, int der, const R *param)    /* K(x)=x^2/c^2 EXP(-x^2/c^2) */
 {
   R c = param[0];
diff --git a/applications/fastsum/kernels.h b/applications/fastsum/kernels.h
index f5766ff7..0742c1b8 100644
--- a/applications/fastsum/kernels.h
+++ b/applications/fastsum/kernels.h
@@ -56,6 +56,7 @@ C kcot(R x, int der, const R *param);                  /**< K(x) = cot(cx) */
 C one_over_cube(R x, int der, const R *param);         /**< K(x) = 1/x^3 */
 C log_sin(R x, int der, const R *param);               /**< K(x) = log(|sin(cx)|) */
 C laplacian_rbf(R x, int der, const R *param);         /**< K(x) = exp(-|x|/c) */
+C der_laplacian_rbf(R x, int der, const R *param);     /**< K(x) = |x|/c exp(-|x|/c) */
 C xx_gaussian(R x, int der, const R *param);           /**< K(x) = x^2/c^2 exp(-x^2/c^2) */
 C absx(R x, int der, const R *param);                  /**< K(x) = |x| */
 /* \} */
diff --git a/julia/fastsum/libfastsumjulia.c b/julia/fastsum/libfastsumjulia.c
index 8b42fd75..c74021b2 100644
--- a/julia/fastsum/libfastsumjulia.c
+++ b/julia/fastsum/libfastsumjulia.c
@@ -55,6 +55,8 @@ int jfastsum_init( fastsum_plan* p, int d, char* s, double* c, unsigned int f, i
 		kernel = log_sin;
 	else if ( strcmp(s, "laplacian_rbf") == 0 )
 		kernel = laplacian_rbf;
+    else if ( strcmp(s, "der_laplacian_rbf") == 0 )
+		kernel = der_laplacian_rbf;
 	else if ( strcmp(s, "xx_gaussian") == 0 )
 		kernel = xx_gaussian;
 	else if ( strcmp(s, "absx") == 0 )
diff --git a/matlab/fastsum/fastsummex.c b/matlab/fastsum/fastsummex.c
index 47fece56..ea50560a 100644
--- a/matlab/fastsum/fastsummex.c
+++ b/matlab/fastsum/fastsummex.c
@@ -134,6 +134,8 @@ static kernel get_kernel(const mxArray *p)
     ker = log_sin;
   else if (strcmp(s, "laplacian_rbf") == 0)
     ker = laplacian_rbf;
+  else if (strcmp(s, "der_laplacian_rbf") == 0)
+    ker = der_laplacian_rbf;
   else if (strcmp(s, "xx_gaussian") == 0)
     ker = xx_gaussian;
   else if (strcmp(s, "absx") == 0)
diff --git a/matlab/fastsum/simple_test.m b/matlab/fastsum/simple_test.m
index d5534f26..6ed15571 100644
--- a/matlab/fastsum/simple_test.m
+++ b/matlab/fastsum/simple_test.m
@@ -39,6 +39,7 @@
 % 'one_over_cube'           K(x) = 1/x^3
 % 'log_sin'                 K(x) = LOG(|SIN(cx)|)
 % 'laplacian_rbf'           K(x) = EXP(-|x|/c)
+% 'der_laplacian_rbf'       K(x) = |x|/c EXP(-|x|/c)
 % 'xx_gaussian'             K(x) = x^2/c^2 EXP(-x^2/c^2) 
 % 'absx'                    K(x) = |x| 
 
diff --git a/matlab/fastsum/test_fastsum.m b/matlab/fastsum/test_fastsum.m
index d316a0c2..ab1aff68 100644
--- a/matlab/fastsum/test_fastsum.m
+++ b/matlab/fastsum/test_fastsum.m
@@ -39,6 +39,7 @@
 % 'one_over_cube'           K(x) = 1/x^3
 % 'log_sin'                 K(x) = LOG(|SIN(cx)|)
 % 'laplacian_rbf'           K(x) = EXP(-|x|/c)
+% 'der_laplacian_rbf'       K(x) = |x|/c EXP(-|x|/c)
 % 'xx_gaussian'             K(x) = x^2/c^2 EXP(-x^2/c^2) 
 % 'absx'                    K(x) = |x|