diff --git a/pub/lib/func.cpp b/pub/lib/func.cpp
index c0ea3bddc43d19e28cb5e9a3e63c58a13f067846..bb5719e0c8cae503b63303eb162239dadd1ec1ca 100644
--- a/pub/lib/func.cpp
+++ b/pub/lib/func.cpp
@@ -71,7 +71,9 @@ double func_exp(double v) { return exp(v); }
 void func_exp(double& r, double& dr, double v, double dv) { r=exp(v); dr=dv*r; }
 
 double func_sin(double v) { return sin(v); }
+void func_sin(double& r, double& dr, double v, double dv) { r=sin(v); dr=dv*fabs(cos(v)); }
 double func_cos(double v) { return cos(v); }
+void func_cos(double& r, double& dr, double v, double dv) { r=cos(v); dr=dv*fabs(sin(v)); }
 double func_tan(double v) { return tan(v); }
 double func_cot(double v) { double t=tan(v); return t==0?0:1/t;}
 
@@ -99,7 +101,7 @@ double func_j0(double v) { return v==0 ? 1 : sin(v)/v; }
 double func_j1(double v) { return v==0 ? 1 : (sin(v)/v-cos(v))/v; }
 
 double func_fac(double v) { return v>-.5 ? gsl_sf_gamma(round(v+1)) : 0; } 
-double func_cata(double v) {
+double func_cata(double v) { // Catalan number ?
     if ( v<=-.5 )
         return 0;
     double x = round(v);
@@ -159,8 +161,7 @@ int func_exit_unless( int a ) { if(!a) exit(0); return 0; }
 //**************************************************************************************************
 
 //double func_pm(double v, double a) { return v; }
-void func_pm(double& r, double& dr, double v, double dv, double a, double da)
-{ r=v; dr=a; }
+void func_pm(double& r, double& dr, double v, double a){ r=v; dr=a; }
 
 template<class A>
 A func_add(A v, A a) { return v+a; }
@@ -569,8 +570,28 @@ void register_fct_e_ee( const char* _tag, func_e_ee _f )
     register_fct_template( _tag, "e", "ee", (void*)_f );
 }
 
+void register_fct_d_ddd( const char* _tag, func_d_ddd _f )
+{
+    register_fct_template( _tag, "d", "ddd", (void*)_f );
+}
+
+void register_fct_e_eee( const char* _tag, func_e_eee _f )
+{
+    register_fct_template( _tag, "e", "eee", (void*)_f );
+}
+
 void NFunctions::initialize(void)
 {
+//  TODO: not yet reimplemented:
+    /*
+    (new CFunc2( "over",func_over, 0,          5   ))->register_me();
+    (new CFunc2( "/",   func_div,  func_div,   6   ))->register_me();
+    (new CFunc2( "mod", func_mod,  0,          6   ))->register_me();
+    (new CFunc2( "div", func_idiv, 0,          6   ))->register_me();
+    (new CFunc2( "xor", func_xor,  0,         10   ))->register_me();
+    (new CFunc3( "?:",  func_cond, func_cond, 13 ))->register_me();
+    */
+    
     // operators by precedence (which should agree with xax_yacc.ypp):
     register_fct_meta ( "+-", 2, "(x,y): x+-y [set error]", 1 );
     register_fct_e_dd ( "+-", func_pm );
@@ -629,187 +650,163 @@ void NFunctions::initialize(void)
     register_fct_i_i  ( "exit_if", func_exit_if );
     register_fct_meta ( "exit_unless", 1, "(x): exit(0) if !x" );
     register_fct_i_i  ( "exit_unless", func_exit_unless );
-    register_fct_meta ( "min", 2, "(x,y): return the smaller of x and y [minimum]" );
-    register_fct_i_ii ( "min", func_min );
-    register_fct_d_dd ( "min", func_min );
-    register_fct_meta ( "max", 2, "(x,y): return the larger of x and y [maximum]" );
-    register_fct_i_ii ( "max", func_max );
-    register_fct_d_dd ( "max", func_max );
-
-//    NFunctions::CastMap.insert( 
-    /*
-    // operators by precedence (as in xax_yacc.ypp):
-    (new CFunc2( "+-",  func_pm,   func_pm,    1   ))->register_me();
-    (new CFunc2( "^",   func_pow,  func_pow,   2   ))->register_me();
-    (new CFunc1( "!",   func_not,  0,          3   ))->register_me();
-    (new CFunc1( "neg", func_neg,  func_neg,   4   ))->register_me();
-    (new CFunc2( "over",func_over, 0,          5   ))->register_me();
-    (new CFunc2( "*",   func_mul,  func_mul,   6   ))->register_me();
-    (new CFunc2( "/",   func_div,  func_div,   6   ))->register_me();
-    (new CFunc2( "mod", func_mod,  0,          6   ))->register_me();
-    (new CFunc2( "div", func_idiv, 0,          6   ))->register_me();
-    (new CFunc2( "+",   func_add,  func_add,   7   ))->register_me();
-    (new CFunc2( "-",   func_sub,  func_sub,   7   ))->register_me();
-    (new CFunc2( ">",   func_gt,   0,          8   ))->register_me();
-    (new CFunc2( ">=",  func_ge,   0,          8   ))->register_me();
-    (new CFunc2( "<",   func_lt,   0,          8   ))->register_me();
-    (new CFunc2( "<=",  func_le,   0,          8   ))->register_me();
-    (new CFunc2( "==",  func_eq,   0,          9   ))->register_me();
-    (new CFunc2( "!=",  func_ne,   0,          9   ))->register_me();
-    (new CFunc2( "xor", func_xor,  0,         10   ))->register_me();
-    (new CFunc2( "&&",  func_and,  0,         11   ))->register_me();
-    (new CFunc2( "||",  func_or ,  0,         12   ))->register_me();
-    (new CFunc3( "?:",  func_cond, func_cond, 13 ))->register_me();
 
     // f( 1 arg )
-    (new CFunc1( "ln",   func_ln,    func_ln, 0,
-           "(x): natural logarithm of x, or 0 if x<=0" ))->register_me();
-    (new CFunc1( "lg",   func_lg,    func_lg, 0,
-           "(x): decadic logarithm of x, or 0 if x<=0" ))->register_me();
-    (new CFunc1( "sqrt", func_sqrt,  func_sqrt, 0,
-           "(x): square root of x, or 0 if x<0" ))->register_me();
-    (new CFunc1( "abs",  func_abs,   func_abs, 0,
-           "(x): absolute value of x" ))->register_me();
-    (new CFunc1( "abs",  func_sign,   func_sign, 0,
-           "(x): sign of x" ))->register_me();
-    (new CFunc1( "exp",  func_exp,   func_exp, 0,
-           "(x): exponential function of x" ))->register_me();
-    (new CFunc1( "sin",  func_sin, 0, 0,
-           "(x): sine of x (where x in radian)" ))->register_me();
-    (new CFunc1( "cos",  func_cos, 0, 0,
-           "(x): cosine of x (where x in radian)" ))->register_me();
-    (new CFunc1( "tan",  func_tan, 0, 0,
-           "(x): tangent of x (where x in radian)" ))->register_me();
-    (new CFunc1( "cot",  func_cot, 0, 0,
-           "(x): cotangent of x (where x in radian)" ))->register_me();
-    (new CFunc1( "sind",  func_sind, 0, 0,
-           "(x): sine of x (where x in degrees)" ))->register_me();
-    (new CFunc1( "cosd",  func_cosd, 0, 0,
-           "(x): cosine of x (where x in degrees)" ))->register_me();
-    (new CFunc1( "tand",  func_tand, 0, 0,
-           "(x): tangent of x (where x in degrees)" ))->register_me();
-    (new CFunc1( "cotd",  func_cotd, 0, 0,
-           "(x): cotangent of x (where x in degrees)" ))->register_me();
-    (new CFunc1( "asin", func_asin, 0, 0,
-           "(x): arc sine of x (result is in radian; return 0 if |x|>1)"
-        ))->register_me();
-    (new CFunc1( "acos", func_acos, 0, 0,
-           "(x): arc cosine of x (result is in radian; return 0 if |x|>1)"
-        ))->register_me();
-    (new CFunc1( "atan", func_atan, 0, 0,
-           "(x): arc tangent of x (result is in radian)" ))->register_me();
-    (new CFunc1( "acot", func_acot, 0, 0,
-           "(x): arc cotangent of x (result is in radian)" ))->register_me();
-    (new CFunc1( "asind", func_asind, 0, 0,
-           "(x): arc sine of x (result is in degrees; return 0 if |x|>1)"
-        ))->register_me();
-    (new CFunc1( "acosd", func_acosd, 0, 0,
-           "(x): arc cosine of x (result is in degrees; return 0 if |x|>1)"
-        ))->register_me();
-    (new CFunc1( "atand", func_atand, 0, 0,
-           "(x): arc tangent of x (result is in degrees)" ))->register_me();
-    (new CFunc1( "acotd", func_acotd, 0, 0,
-           "(x): arc cotangent of x (result is in degrees)" ))->register_me();
-    (new CFunc1( "sinh", func_sinh, 0, 0,
-           "(x): hyperbolic sine of x" ))->register_me();
-    (new CFunc1( "cosh", func_cosh, 0, 0,
-           "(x): hyperbolic cosine of x" ))->register_me();
-    (new CFunc1( "tanh", func_tanh, 0, 0,
-           "(x): hyperbolic tangent of x" ))->register_me();
-    (new CFunc1( "coth", func_coth, 0, 0,
-           "(x): hyperbolic tangent of x" ))->register_me();
-    (new CFunc1( "j0", func_j0, 0, 0,
-           "(x): spherical Bessel function x" ))->register_me();
-    (new CFunc1( "j1", func_j1, 0, 0,
-           "(x): spherical Bessel function x" ))->register_me();
-    (new CFunc1( "fac",func_fac, 0, 0,
-           "(x): factorial of nearest integer of x" ))->register_me();
-    (new CFunc1( "cata",func_cata, 0, 0,
-           "(x): Catalan number of nearest integer of x" ))->register_me();
-    (new CFunc1( "gamma",func_gamma, 0, 0,
-           "(x): gamma function of x (i.e. factorial of x-1)" ))->register_me();
-    (new CFunc1( "lgamma",func_lgamma, 0, 0,
-           "(x): ln of gamma of x" ))->register_me();
-    (new CFunc1( "debye1",func_debye1, 0, 0,
-           "(x): Debye function D1 of x" ))->register_me();
-    (new CFunc1( "debye3",func_debye3, 0, 0,
-           "(x): Debye function D3 of x" ))->register_me();
-    (new CFunc1( "debyeu2",func_debye_msd, 0, 0,
-           "(x): Debye mean squared displacement's temperature dependence"
-        ))->register_me();
-    (new CFunc1( "debyeui",func_debye_ui, 0, 0,
-           "(x): Debye internal energy" ))->register_me();
-    (new CFunc1( "debyecv",func_debye_cv, 0, 0,
-           "(x): Debye heat capacity" ))->register_me();
-    (new CFunc1( "erfP", func_erfP, 0, 0,
-           "(x): ?" ))->register_me();
-    (new CFunc1( "erfQ", func_erfQ , 0, 0,
-           "(x): ?"  ))->register_me();
-    (new CFunc1( "erf",  func_erf, 0, 0,
-           "(x): error function of x" ))->register_me();
-    (new CFunc1( "erfc", func_erfc, 0, 0,
-           "(x): complementary error function of x" ))->register_me();
-    (new CFunc1( "dawson", func_dawson, 0, 0,
-           "(x): Dawson function of x" ))->register_me();
-    (new CFunc1( "sinc", func_sinc, 0, 0,
-           "(x): sinus cardinalis, sin(x)/x" ))->register_me();
-    (new CFunc1( "ceil", func_ceil, 0, 0,
-           "(x): the smallest integer i with i>=x" ))->register_me();
-    (new CFunc1( "floor",func_floor, 0, 0,
-           "(x): the largest integer i with i<=x" ))->register_me();
-    (new CFunc1( "nint", func_nint, 0, 0,
-           "(x): the integer nearest to x" ))->register_me();
-
-    (new CFunc1( "diehl", func_diehl, 0, 0,
-           "(cauchywid/gausswid): normalized convolution gauss(*)cauchy"        ))->register_me();
-    (new CFunc1( "lndiehl", func_lndiehl, 0, 0,
-           "(cauchywid/gausswid): ln of normalized convolution gauss(*)cauchy"        ))->register_me();
-
-    (new CFunc1( "exit", func_exit, 0, 0,
-           "(x): exit with return value nint(x)" ))->register_me();
-
+    register_fct_meta ( "ln", 1, "(x): natural logarithm of x, or 0 if x<=0" );
+    register_fct_d_d  ( "ln", func_ln );
+    register_fct_e_e  ( "ln", func_ln );
+    register_fct_meta ( "lg", 1, "(x): decadic logarithm of x, or 0 if x<=0" );
+    register_fct_d_d  ( "lg", func_lg );
+    register_fct_e_e  ( "lg", func_lg );
+    register_fct_meta ( "sqrt", 1, "(x): square root of x, or 0 if x<0" );
+    register_fct_d_d  ( "sqrt", func_sqrt );
+    register_fct_e_e  ( "sqrt", func_sqrt );
+    register_fct_meta ( "abs", 1, "(x): absolute value of x" );
+    register_fct_i_i  ( "abs", func_abs );
+    register_fct_d_d  ( "abs", func_abs );
+    register_fct_e_e  ( "abs", func_abs );
+    register_fct_meta ( "sign", 1, "(x): sign of x" );
+    register_fct_i_i  ( "sign", func_sign );
+    register_fct_d_d  ( "sign", func_sign );
+    register_fct_e_e  ( "sign", func_sign );
+    register_fct_meta ( "exp", 1, "(x): exponential function of x" );
+    register_fct_d_d  ( "exp", func_exp );
+    register_fct_e_e  ( "exp", func_exp );
+    register_fct_meta ( "sin", 1, "(x): sine of x (where x in radian)" );
+    register_fct_d_d  ( "sin", func_sin );
+    register_fct_e_e  ( "sin", func_sin );
+    register_fct_meta ( "cos", 1, "(x): cosine of x (where x in radian)" );
+    register_fct_d_d  ( "cos", func_cos );
+    register_fct_e_e  ( "cos", func_cos );
+    register_fct_meta ( "tan", 1, "(x): tangent of x (where x in radian)" );
+    register_fct_d_d  ( "tan", func_tan );
+    register_fct_meta ( "cot", 1, "(x): cotangent of x (where x in radian)" );
+    register_fct_d_d  ( "cot", func_cot );
+    register_fct_meta ( "sind", 1, "(x): sine of x (where x in degrees)" );
+    register_fct_d_d  ( "sind", func_sind );
+    register_fct_meta ( "cosd", 1, "(x): cosine of x (where x in degrees)" );
+    register_fct_d_d  ( "cosd", func_cosd );
+    register_fct_meta ( "tand", 1, "(x): tangent of x (where x in degrees)" );
+    register_fct_d_d  ( "tand", func_tand );
+    register_fct_meta ( "cotd", 1, "(x): cotangent of x (where x in degrees)" );
+    register_fct_d_d  ( "cotd", func_cotd );
+    register_fct_meta ( "asin", 1, "(x): arc sine of x (result is in radian; 0 if |x|>1)" );
+    register_fct_d_d  ( "asin", func_asin );
+    register_fct_meta ( "acos", 1, "(x): arc cosine of x (result is in radian; 0 if |x|>1)" );
+    register_fct_d_d  ( "acos", func_acos );
+    register_fct_meta ( "atan", 1, "(x): arc tangent of x (result is in radian)" );
+    register_fct_d_d  ( "atan", func_atan );
+    register_fct_meta ( "acot", 1, "(x): arc cotangent of x (result is in radian)" );
+    register_fct_d_d  ( "acot", func_acot );
+    register_fct_meta ( "asind", 1, "(x): arc sine of x (result is in degrees; 0 if |x|>1)" );
+    register_fct_d_d  ( "asind", func_asind );
+    register_fct_meta ( "acosd", 1, "(x): arc cosine of x (result is in degrees; 0 if |x|>1)" );
+    register_fct_d_d  ( "acosd", func_acosd );
+    register_fct_meta ( "atand", 1, "(x): arc tangent of x (result is in degrees)" );
+    register_fct_d_d  ( "atand", func_atand );
+    register_fct_meta ( "acotd", 1, "(x): arc cotangent of x (result is in degrees)" );
+    register_fct_d_d  ( "acotd", func_acotd );
+    register_fct_meta ( "sinh", 1, "(x): hyperbolic sine of x" );
+    register_fct_d_d  ( "sinh", func_sinh );
+    register_fct_meta ( "cosh", 1, "(x): hyperbolic cosine of x" );
+    register_fct_d_d  ( "cosh", func_cosh );
+    register_fct_meta ( "tanh", 1, "(x): hyperbolic tangent of x" );
+    register_fct_d_d  ( "tanh", func_tanh );
+    register_fct_meta ( "coth", 1, "(x): hyperbolic tangent of x" );
+    register_fct_d_d  ( "coth", func_coth );
+    register_fct_meta ( "j0", 1, "(x): spherical Bessel function x" );
+    register_fct_d_d  ( "j0", func_j0 );
+    register_fct_meta ( "j1", 1, "(x): spherical Bessel function x" );
+    register_fct_d_d  ( "j1", func_j1 );
+    register_fct_meta ( "fac", 1, "(x): factorial of nearest integer of x" );
+    register_fct_d_d  ( "fac", func_fac );
+    register_fct_meta ( "cata", 1, "(x): Catalan number of nearest integer of x" );
+    register_fct_d_d  ( "cata", func_cata );
+    register_fct_meta ( "gamma", 1, "(x): gamma function of x (i.e. factorial of x-1)" );
+    register_fct_d_d  ( "gamma", func_gamma );
+    register_fct_meta ( "lgamma", 1, "(x): ln of gamma of x" );
+    register_fct_d_d  ( "lgamma", func_lgamma );
+    register_fct_meta ( "debye1", 1, "(x): Debye function D1 of x" );
+    register_fct_d_d  ( "debye1", func_debye1 );
+    register_fct_meta ( "debye3", 1, "(x): Debye function D3 of x" );
+    register_fct_d_d  ( "debye3", func_debye3 );
+    register_fct_meta ( "debyeu2", 1, "(x): Debye mean squared displacement's temperature dependence" );
+    register_fct_d_d  ( "debyeu2", func_debye_msd );
+    register_fct_meta ( "debyeui", 1, "(x): Debye internal energy" );
+    register_fct_d_d  ( "debyeui", func_debye_ui );
+    register_fct_meta ( "debyecv", 1, "(x): Debye heat capacity" );
+    register_fct_d_d  ( "debyecv", func_debye_cv );
+    register_fct_meta ( "erfP", 1, "(x): ?" );
+    register_fct_d_d  ( "erfP", func_erfP );
+    register_fct_meta ( "erfQ", 1, "(x): ?" );
+    register_fct_d_d  ( "erfQ", func_erfQ );
+    register_fct_meta ( "erf", 1, "(x): error function of x" );
+    register_fct_d_d  ( "erf", func_erf );
+    register_fct_meta ( "erfc", 1, "(x): complementary error function of x" );
+    register_fct_d_d  ( "erfc", func_erfc );
+    register_fct_meta ( "dawson", 1, "(x): Dawson function of x" );
+    register_fct_d_d  ( "dawson", func_dawson );
+    register_fct_meta ( "sinc", 1, "(x): sinus cardinalis, sin(x)/x" );
+    register_fct_d_d  ( "sinc", func_sinc );
+    register_fct_meta ( "ceil", 1, "(x): the smallest integer i with i>=x" );
+    register_fct_d_d  ( "ceil", func_ceil );
+    register_fct_meta ( "floor", 1, "(x): the largest integer i with i<=x" );
+    register_fct_d_d  ( "floor", func_floor );
+    register_fct_meta ( "nint", 1, "(x): the integer nearest to x" );
+    register_fct_d_d  ( "nint", func_nint );
+
+    register_fct_meta ( "diehl", 1, "(cauchywid/gausswid): normalized convolution gauss(*)cauchy" );
+    register_fct_d_d  ( "diehl", func_diehl );
+    register_fct_meta ( "lndiehl", 1, "(cauchywid/gausswid): ln of normalized convolution gauss(*)cauchy" );
+    register_fct_d_d  ( "lndiehl", func_lndiehl );
+
+    register_fct_meta ( "exit", 1, "(x): exit with return value nint(x)" );
+    register_fct_i_i  ( "exit", func_exit );
+    register_fct_meta ( "exit_if", 1, "(x): exit(0) if x" );
+    register_fct_i_i  ( "exit_if", func_exit_if );
+    register_fct_meta ( "exit_unless", 1, "(x): exit(0) if !x" );
+    register_fct_i_i  ( "exit_unless", func_exit_unless );
 
     // f(2 args)
-    (new CFunc2( "min2",  func_min, 0, 0,
-           "(x,y): the smaller of the two arguments x and y" ))->register_me();
-    (new CFunc2( "max2",  func_max, 0, 0,
-           "(x,y): the smaller of the two arguments x and y" ))->register_me();
-    (new CFunc2( "ran",  func_ran, 0, 0,
-           "(x,y): a random number between x and y" ))->register_me();
-
-    (new CFunc2( "gauss",   func_gauss, 0, 0,
-           "(x,s): the normalized Gaussian exp(-x^2/2/s^2)/sqrt(2 pi)/s"
-        ))->register_me();
-    (new CFunc2( "gnn",     func_gaussnn, 0, 0,
-           "(x,s): the unnormalized Gaussian exp(-x^2/2/s^2)" ))->register_me();
-    (new CFunc2( "cauchy",  func_cauchy, 0, 0,
-           "(x,w): the Cauchy-Lorentz function ?/(x^2+w^2)" ))->register_me();
-
+    register_fct_meta ( "min2", 2, "(x,y): the smaller of the two arguments x and y" );
+    register_fct_i_ii ( "min2", func_min );
+    register_fct_d_dd ( "min2", func_min );
+    register_fct_meta ( "max2", 2, "(x,y): the smaller of the two arguments x and y" );
+    register_fct_i_ii ( "max2", func_max );
+    register_fct_d_dd ( "max2", func_max );
+    register_fct_meta ( "ran", 2, "(x,y): a random number between x and y" );
+    register_fct_d_dd ( "ran", func_ran );
+
+    register_fct_meta ( "gauss", 2, "(x,s): the normalized Gaussian exp(-x^2/2/s^2)/sqrt(2 pi)/s" );
+    register_fct_d_dd ( "gauss", func_gauss );
+    register_fct_meta ( "gnn", 2, "(x,s): the unnormalized Gaussian exp(-x^2/2/s^2)" );
+    register_fct_d_dd ( "gnn", func_gaussnn );
+    register_fct_meta ( "cauchy", 2, "(x,w): the Cauchy-Lorentz function ?/(x^2+w^2)" );
+    register_fct_d_dd ( "cauchy", func_cauchy );
 
     // f(3 args)
-    (new CFunc3( "rehavneg", func_re_havneg, 0, 0,
-           "(x,y,z): real part of the Havriliak-Negami function"        ))->register_me();
-    (new CFunc3( "imhavneg", func_im_havneg, 0, 0,
-           "(x,y,z): imaginary part of the Havriliak-Negami function"        ))->register_me();
-    (new CFunc3( "q4w",      func_q4w, 0, 0,
-           "(x,y,z): ?" ))->register_me();
-    (new CFunc3( "cauchy2",  func_cauchy2, 0, 0,
-           "(x,y,z): ?" ))->register_me();
-    (new CFunc3( "kwwc", func_kwwc, 0, 0,
-           "(w,tau,b): Fourier cosine transform (t->w) of exp((t/tau)^b)"        ))->register_me();
-    (new CFunc3( "kwws", func_kwws, 0, 0,
-           "(w,tau,b): Fourier sine transform (t->w) of exp((t/tau)^b)"        ))->register_me();
-    (new CFunc3( "kwwp", func_kwwp, 0, 0,
-           "(w,tau,b): primitive of Fourier cosine transform (t->w) of exp((t/tau)^b)"        ))->register_me();
-    (new CFunc3( "voigt", func_voigt, 0, 0,
-           "(x,sigma,gamma): convolution of Gaussian(x,sigma) and Lorentzian(x,gamma)"        ))->register_me();
-    (new CFunc3( "zorn", func_zorn, 0, 0,
-           "(I,<I>,s): Zorn's multiple-scattering corrected elastic intensity"        ))->register_me();
-    (new CFunc3( "zorn2", func_zorn_gauss, 0, 0,
-           "(q,<u^2>,s): Zorn's multiple-scattering corrected Gaussian elastic intensity for Si111" ))->register_me();
-    (new CFunc3( "rrdm", func_rrdm, 0, 0,
-           "(w*t0,EA_mean/T,EA_stdv/EA_mean: rotational rate distribution model"        ))->register_me();
-    (new CFunc3( "rotdiff", func_rotdiff, 0, 0,
-           "(w,tau,qb: rotational diffusion spectrum)"))->register_me();
-    */
+    register_fct_meta ( "rehavneg", 3, "(x,y,z): real part of the Havriliak-Negami function" );
+    register_fct_d_ddd( "rehavneg", func_re_havneg );
+    register_fct_meta ( "imhavneg", 3, "(x,y,z): imaginary part of the Havriliak-Negami function" );
+    register_fct_d_ddd( "imhavneg", func_im_havneg );
+    register_fct_meta ( "q4w", 3, "(x,y,z): ?" );
+    register_fct_d_ddd( "q4w", func_q4w );
+    register_fct_meta ( "cauchy2", 3, "(x,y,z): ?" );
+    register_fct_d_ddd( "cauchy2", func_cauchy2 );
+    register_fct_meta ( "kwwc", 3, "(w,tau,b): Fourier cosine transform (t->w) of exp((t/tau)^b)" );
+    register_fct_d_ddd( "kwwc", func_kwwc );
+    register_fct_meta ( "kwws", 3, "(w,tau,b): Fourier sine transform (t->w) of exp((t/tau)^b)" );
+    register_fct_d_ddd( "kwws", func_kwws );
+    register_fct_meta ( "kwwp", 3, "(w,tau,b): primitive of Fourier cosine transform (t->w) of exp((t/tau)^b)" );
+    register_fct_d_ddd( "kwwp", func_kwwp );
+    register_fct_meta ( "voigt", 3, "(x,sigma,gamma): convolution of Gaussian(x,sigma) and Lorentzian(x,gamma)" );
+    register_fct_d_ddd( "voigt", func_voigt );
+    register_fct_meta ( "zorn", 3, "(I,<I>,s): Zorn's multiple-scattering corrected elastic intensity" );
+    register_fct_d_ddd( "zorn", func_zorn );
+    register_fct_meta ( "zorn2", 3, "(q,<u^2>,s): Zorn's multiple-scattering corrected Gaussian elastic intensity for Si111" );
+    register_fct_d_ddd( "zorn2", func_zorn_gauss );
+    register_fct_meta ( "rrdm", 3, "(w*t0,EA_mean/T,EA_stdv/EA_mean: rotational rate distribution model" );
+    register_fct_d_ddd( "rrdm", func_rrdm );
+    register_fct_meta ( "rotdiff", 3, "(w,tau,qb: rotational diffusion spectrum)" );
+    register_fct_d_ddd( "rotdiff", func_rotdiff );
 }
diff --git a/pub/lib/func.hpp b/pub/lib/func.hpp
index 40ca96f5717adfd1e046ce82aa9d188b53066e32..df59fac833602651712237f09f7bf412c891a75d 100644
--- a/pub/lib/func.hpp
+++ b/pub/lib/func.hpp
@@ -23,22 +23,12 @@ typedef double (*func_d_ii) (int, int);
 typedef double (*func_d_id) (int, double);
 typedef double (*func_d_di) (double, int);
 typedef double (*func_d_dd) (double, double);
-typedef void   (*func_e_dd) (double&,double&, double,double, double,double);
+typedef void   (*func_e_dd) (double&,double&, double, double);
 typedef void   (*func_e_ee) (double&,double&, double,double, double,double);
 
+typedef double (*func_d_ddd) (double, double, double);
 typedef void   (*func_e_eee) (double&,double&, double,double, double,double, double,double);
 
-typedef double (*func_f1) (double);
-typedef double (*func_f2) (double, double);
-typedef double (*func_f3) (double, double, double);
-
-//  Derivatives thereof.
-
-typedef void (*deri_f1) (double&,double&,double,double);
-typedef void (*deri_f2) (double&,double&,double,double,double,double);
-typedef void (*deri_f3) (double&,double&,
-                         double,double,double,double,double,double);
-
 //! A wrapper holding a function with given input and output types.
 
 class CTypedFunc {
diff --git a/pub/lib/node.cpp b/pub/lib/node.cpp
index 21e2b9ea4e8e59185bf2aa056200ea571157e888..d60e0148b625687f0a92bb224f0764fb67841bc3 100644
--- a/pub/lib/node.cpp
+++ b/pub/lib/node.cpp
@@ -62,18 +62,18 @@ const RObj CNodeFun::tree_val( const CContext& ctx ) const
     try {
          // evaluate arguments, determine if all results are scalar.
         vector<RObj> pa(narg);
-        string base_types = "";
+        string basetypes = "";
         bool is_scalar = true;
         for ( int iarg=0; iarg<narg; ++iarg ) {
             pa[iarg] = arg[iarg]->tree_val( ctx );
-            base_types += pa[iarg]->base_type();
+            basetypes += pa[iarg]->base_type();
             if ( pa[iarg]->vectorial )
                 is_scalar = false;
         }
         
-        const CTypedFunc* tf = fu->find_tyfu( base_types, ctx.want_error );
+        const CTypedFunc* tf = fu->find_tyfu( basetypes, ctx.want_error );
         if( !tf )
-            throw "function " + fu->tag + " not implemented for argument types " + base_types;
+            throw "function " + fu->tag + " not implemented for argument types " + basetypes;
         void* f = tf->f;
 
         if( is_scalar ) {
@@ -88,7 +88,7 @@ const RObj CNodeFun::tree_val( const CContext& ctx ) const
                 else if  ( tf->intypes=="dd" )
                     ret = (*(func_i_dd)(f))( pa[0]->to_r(), pa[1]->to_r() );
                 else
-                    throw "BUG: unexpected intypes";
+                    throw "BUG: unexpected intypes '"+basetypes+"', cast to '"+tf->intypes+"'";
                 return PObjInt( new CObjInt( ret ) );
 
             } else if ( tf->outtype=="d" ) {
@@ -101,14 +101,18 @@ const RObj CNodeFun::tree_val( const CContext& ctx ) const
                     ret = (*(func_d_ii)(f))( pa[0]->to_i(), pa[1]->to_i() );
                 else if  ( tf->intypes=="dd" )
                     ret = (*(func_d_dd)(f))( pa[0]->to_r(), pa[1]->to_r() );
+                else if  ( tf->intypes=="ddd" )
+                    ret = (*(func_d_ddd)(f))( pa[0]->to_r(), pa[1]->to_r(), pa[2]->to_r() );
                 else
-                    throw "BUG: unexpected intypes";
+                    throw "BUG: unexpected intypes '"+basetypes+"', cast to '"+tf->intypes+"'";
                 return PObjDble( new CObjDble( ret ) );
             } else if ( tf->outtype=="e" ) {
                 double val, err;
                 if       ( tf->intypes=="e" )
                     (*(func_e_e)(f))( val, err,
                              pa[0]->to_r(), pa[0]->to_dr() );
+                else if  ( tf->intypes=="dd" )
+                    (*(func_e_dd)(f))( val, err, pa[0]->to_r(), pa[1]->to_r() );
                 else if  ( tf->intypes=="ee" )
                     (*(func_e_ee)(f))( val, err,
                              pa[0]->to_r(), pa[0]->to_dr(),
@@ -119,10 +123,10 @@ const RObj CNodeFun::tree_val( const CContext& ctx ) const
                              pa[1]->to_r(), pa[1]->to_dr(),
                              pa[2]->to_r(), pa[2]->to_dr() );
                 else
-                    throw "BUG: unexpected intypes";
+                    throw "BUG: unexpected intypes '"+basetypes+"', cast to '"+tf->intypes+"'";
                 return PObjDble( new CObjDble( val, err ) );
             } else
-                throw "BUG: unexpected outtype";
+                throw "BUG: unexpected outtype '"+tf->outtype+"' after cast";
         }
         
         throw "incomplete";
@@ -137,15 +141,6 @@ const RObj CNodeFun::tree_val( const CContext& ctx ) const
             else if ( pan[iarg] = PCAST<const CObjNumber>(pa[iarg]) )
  
         if( is_scalar ) {
-            PObjDble pret( new CObjDble() );
-            if ( want_err )
-                (*(fu->d3))( pret->val, pret->err,
-                             pan[0]->to_r(), pan[0]->to_dr(),
-                             pan[1]->to_r(), pan[1]->to_dr(),
-                             pan[2]->to_r(), pan[2]->to_dr() );
-            else
-                pret->val = (*(fu->f3))( pan[0]->to_r(), pan[1]->to_r(), pan[2]->to_r() );
-            return pret;
         } else {
             int n=0;
             for ( int iarg=0; iarg<narg; ++iarg ) {
diff --git a/pub/lib/obj.hpp b/pub/lib/obj.hpp
index 87cd427641b2d9b1b452a7846080ed9ad3ea5cff..7dc5d3d8871500772d090ba71e164050da969bf5 100644
--- a/pub/lib/obj.hpp
+++ b/pub/lib/obj.hpp
@@ -59,7 +59,7 @@ class CObjDble : public CObjNumber {
  public:
     double val;         //!< Scalar value.
     double err;         //!< Error of scalar value.
-    CObjDble( double _val, double _err=0 ) : CObjNumber(), val(_val), err(_err) {};
+    CObjDble( double _val, double _err=-INFINITY ) : CObjNumber(), val(_val), err(_err) {};
     bool has_err() const { return err>=0; };
     char base_type() const { return has_err() ? 'e' : 'd'; };
     string result_info() const;