Fixes #728

src/Classic/AbsBeamline/EndFieldModel/AsymmetricEnge.h

@@ -142,8 +142,19 @@ void AsymmetricEnge::setX0End(double x0) {
 }
 
 double AsymmetricEnge::function(double x, int n) const {
-    return engeStart_m->getEnge(x-engeStart_m->getX0(), n)+
-           engeEnd_m->getEnge(-x-engeEnd_m->getX0(), n)-1;
+    // f(x) = E(x-x0) + E(-x-x0) - 1
+    // f^{(2n)} = E^{(2n)}(x-x0) + E^{(2n)}(-x-x0)
+    // f^{(2n+1)} = E^{(2n)}(x-x0) - E^{(2n)}(-x-x0)
+    if (n == 0) {
+        return engeStart_m->getEnge(x-engeStart_m->getX0(), n)+
+               engeEnd_m->getEnge(-x-engeEnd_m->getX0(), n)-1;
+    } else if (n%2) {
+        return engeStart_m->getEnge(x-engeStart_m->getX0(), n)-
+               engeEnd_m->getEnge(-x-engeEnd_m->getX0(), n);
+    } else {
+        return engeStart_m->getEnge(x-engeStart_m->getX0(), n)+
+               engeEnd_m->getEnge(-x-engeEnd_m->getX0(), n);
+    }
 }
 
 AsymmetricEnge* AsymmetricEnge::clone() const {

tests/classic_src/AbsBeamline/EndFieldModel/AsymmetricEngeTest.cpp

@@ -25,3 +25,20 @@ TEST(AsymmetricEngeTest, FunctionTest) {
     }
 
 }
+
+TEST(AsymmetricEngeTest, DerivativeTest) {
+    endfieldmodel::AsymmetricEnge enge({0.0, 1.0}, 10.0, 3.0, 
+                                       {0.0, 4.0}, 11.0, 6.0);
+    double ds = 0.01;
+    for (int n = 1; n < 5; ++n) {
+        for (double s = -20.0; s < 20.0; s += 1.0) {
+            double fplus = enge.function(s+ds, n-1);
+            double fminus = enge.function(s-ds, n-1);
+            double test = enge.function(s, n);
+            double ref = (fplus-fminus)/2/ds;
+            double delta = abs(ref*1e-3);
+            EXPECT_NEAR(test, ref, delta) << s << " " << test << " " << ref;
+        }
+    }
+}
+