diff --git a/Project.toml b/Project.toml
index ec638c95..830f72a5 100644
--- a/Project.toml
+++ b/Project.toml
@@ -1,7 +1,7 @@
 name = "SummationByPartsOperators"
 uuid = "9f78cca6-572e-554e-b819-917d2f1cf240"
 author = ["Hendrik Ranocha"]
-version = "0.5.43"
+version = "0.5.44"
 
 [deps]
 ArgCheck = "dce04be8-c92d-5529-be00-80e4d2c0e197"
diff --git a/dev/Mattsson2014_dev.jl b/dev/Mattsson2014_dev.jl
index 28de36d2..097f3af5 100644
--- a/dev/Mattsson2014_dev.jl
+++ b/dev/Mattsson2014_dev.jl
@@ -36,7 +36,7 @@ R = vcat(r1, r2, r3, r4, r5, r6, r7)
 n1 = [13//10 -12//5 9//10 1//5 0 0 0]
 n2 = [n1[2] 26//5 -16//5 2//5 0 0 0]
 n3 = [n1[3] n2[3] 47//10 -17//5 1 0 0]
-n4 = [n1[4] n2[4] n3[4] -29//5 -4 1 0]
+n4 = [n1[4] n2[4] n3[4] 29//5 -4 1 0]
 n5 = [0 0 1//1 -4 6 -4 1]
 n6 = [0 0 0 1//1 -4 6 -4]
 n7 = [0 0 0 0 1//1 -4 6]
diff --git a/src/SBP_coefficients/Mattsson2014.jl b/src/SBP_coefficients/Mattsson2014.jl
index e456028b..2faf19db 100644
--- a/src/SBP_coefficients/Mattsson2014.jl
+++ b/src/SBP_coefficients/Mattsson2014.jl
@@ -4,7 +4,7 @@
 
 Coefficients of the SBP operators given in
 - Mattsson (2014)
-  Diagonal-norm summation by parts operators for fiite difference approximations
+  Diagonal-norm summation by parts operators for finite difference approximations
     of third and fourth derivatives.
   Journal of Computational Physics 274, pp. 432-454.
 """
@@ -16,7 +16,7 @@ function Base.show(io::IO, source::Mattsson2014)
     else
         print(io,
             "Mattsson (2014) \n",
-            "  Diagonal-norm summation by parts operators for fiite difference approximations\n",
+            "  Diagonal-norm summation by parts operators for finite difference approximations\n",
             "    of third and fourth derivatives. \n",
             "  Journal of Computational Physics 274, pp. 432-454.")
     end
@@ -458,7 +458,8 @@ function third_derivative_coefficients(source::Mattsson2014, order::Int, T=Float
         upper_coef = SVector(T(-1), T(1//2))
         central_coef = T(0)
         lower_coef = -upper_coef
-        left_weights = SVector(T(1//2))
+        # left_weights = SVector(T(1//2))
+        left_weights = SVector(T(1//2), T(1), T(1))
         right_weights = left_weights
         left_boundary_derivatives = (
             # first derivative
@@ -712,7 +713,7 @@ function fourth_derivative_coefficients(source::Mattsson2014, order::Int, T=Floa
             DerivativeCoefficientRow{T,1,6}(SVector(T(1//5),
                                                     T(2//5),
                                                     T(-17//5),
-                                                    T(-29//5),
+                                                    T(29//5),
                                                     T(-4),
                                                     T(1) )),
         )
@@ -720,7 +721,7 @@ function fourth_derivative_coefficients(source::Mattsson2014, order::Int, T=Floa
         upper_coef = SVector(T(-4), T(1))
         central_coef = T(6)
         lower_coef = upper_coef
-        left_weights = SVector(T(1//2))
+        left_weights = SVector(T(1//2), T(1), T(1), T(1))
         right_weights = left_weights
         left_boundary_derivatives = (
             # first derivative
diff --git a/test/SBP_operators_test.jl b/test/SBP_operators_test.jl
index 4441585d..fb17e345 100644
--- a/test/SBP_operators_test.jl
+++ b/test/SBP_operators_test.jl
@@ -467,7 +467,7 @@ for source in D_test_list, T in (Float32,Float64)
 end
 
 # Accuracy tests of third derivative operators.
-for source in D_test_list, T in (Float32,Float64)
+@testset "third-derivative operators" for source in D_test_list, T in (Float32,Float64)
     xmin = -one(T)
     xmax = 2*one(T)
     N = 101
@@ -499,16 +499,12 @@ for source in D_test_list, T in (Float32,Float64)
         mul!(res, D, x0)
         @test all(i->abs(res[i]) < eps(T), eachindex(res))
         mul!(res, D, x1)
-        @test all(i->abs(res[i]) < 5000*eps(T), eachindex(res))
-        mul!(res, D, x2)
-        @test all(i->abs(res[i]) < 5000*eps(T), eachindex(res))
+        @test all(i->abs(res[i]) < 50_000*eps(T), eachindex(res))
         # only interior
         mul!(res, D, x2)
-        @test all(i->abs(res[i]) < 5000*eps(T), inner_indices)
+        @test all(i->abs(res[i]) < 160_000*eps(T), inner_indices)
         mul!(res, D, x3)
-        @test all(i->res[i] ≈ 6*x0[i], inner_indices)
-        mul!(res, D, x4)
-        @test any(i->!(res[i] ≈ 24*x1[i]), inner_indices)
+        @test all(i->abs(res[i] - 6*x0[i]) < 340_000*eps(T), inner_indices)
         # boundary: first derivative
         @test abs(derivative_left( D, x0, Val{1}())) < eps(T)
         @test abs(derivative_right(D, x0, Val{1}())) < eps(T)
@@ -519,8 +515,8 @@ for source in D_test_list, T in (Float32,Float64)
         # boundary: second derivative
         @test abs(derivative_left( D, x0, Val{2}())) < eps(T)
         @test abs(derivative_right(D, x0, Val{2}())) < eps(T)
-        @test abs(derivative_left( D, x1, Val{2}())) < eps(T)
-        @test abs(derivative_right(D, x1, Val{2}())) < eps(T)
+        @test abs(derivative_left( D, x1, Val{2}())) < 600 * eps(T)
+        @test abs(derivative_right(D, x1, Val{2}())) < 600 * eps(T)
         @test derivative_left( D, x2, Val{2}()) ≈ 2
         @test derivative_right(D, x2, Val{2}()) ≈ 2
     end
@@ -648,7 +644,7 @@ for source in D_test_list, T in (Float32,Float64)
 end
 
 # Accuracy tests of fourth derivative operators.
-for source in D_test_list, T in (Float32,Float64)
+@testset "fourth-derivative operators" for source in D_test_list, T in (Float32,Float64)
     xmin = -one(T)
     xmax = 2*one(T)
     N = 101
@@ -678,18 +674,16 @@ for source in D_test_list, T in (Float32,Float64)
         @test issymmetric(D) == false
         # interior and boundary
         mul!(res, D, x0)
-        @test all(i->abs(res[i]) < eps(T), eachindex(res))
+        @test all(i->abs(res[i]) < 10 * eps(T) / D.Δx^4, eachindex(res))
         mul!(res, D, x1)
-        @test all(i->abs(res[i]) < 5000*eps(T), eachindex(res))
+        @test all(i->abs(res[i]) < 100 * eps(T) / D.Δx^4, eachindex(res))
         mul!(res, D, x2)
-        @test all(i->abs(res[i]) < 5000*eps(T), eachindex(res))
+        @test all(i->abs(res[i]) < 100 * eps(T) / D.Δx^4, eachindex(res))
         # only interior
         mul!(res, D, x2)
-        @test all(i->abs(res[i]) < 5000*eps(T), inner_indices)
+        @test all(i->abs(res[i]) < 100 * eps(T) / D.Δx^4, inner_indices)
         mul!(res, D, x3)
-        @test all(i->abs(res[i]) < 5000*eps(T), inner_indices)
-        mul!(res, D, x4)
-        @test any(i->!(res[i] ≈ 24*x0[i]), inner_indices)
+        @test all(i->abs(res[i]) < 100 * eps(T) / D.Δx^4, inner_indices)
         # boundary: first derivative
         @test abs(derivative_left( D, x0, Val{1}())) < eps(T)
         @test abs(derivative_right(D, x0, Val{1}())) < eps(T)
@@ -700,17 +694,17 @@ for source in D_test_list, T in (Float32,Float64)
         # boundary: second derivative
         @test abs(derivative_left( D, x0, Val{2}())) < eps(T)
         @test abs(derivative_right(D, x0, Val{2}())) < eps(T)
-        @test abs(derivative_left( D, x1, Val{2}())) < eps(T)
-        @test abs(derivative_right(D, x1, Val{2}())) < eps(T)
-        @test derivative_left( D, x2, Val{2}()) ≈ one(T)
-        @test derivative_right(D, x2, Val{2}()) ≈ one(T)
+        @test abs(derivative_left( D, x1, Val{2}())) < 600 * eps(T)
+        @test abs(derivative_right(D, x1, Val{2}())) < 600 * eps(T)
+        @test derivative_left( D, x2, Val{2}()) ≈ 2 * one(T)
+        @test derivative_right(D, x2, Val{2}()) ≈ 2 * one(T)
         # boundary: third derivative
-        @test abs(derivative_left( D, x0, Val{3}())) < eps(T)
-        @test abs(derivative_right(D, x0, Val{3}())) < eps(T)
-        @test abs(derivative_left( D, x1, Val{3}())) < eps(T)
-        @test abs(derivative_right(D, x1, Val{3}())) < eps(T)
-        @test abs(derivative_left( D, x2, Val{3}())) < eps(T)
-        @test abs(derivative_right(D, x2, Val{3}())) < eps(T)
+        @test abs(derivative_left( D, x0, Val{3}())) < 10 * eps(T) / D.Δx^3
+        @test abs(derivative_right(D, x0, Val{3}())) < 10 * eps(T) / D.Δx^3
+        @test abs(derivative_left( D, x1, Val{3}())) < 10 * eps(T) / D.Δx^3
+        @test abs(derivative_right(D, x1, Val{3}())) < 10 * eps(T) / D.Δx^3
+        @test abs(derivative_left( D, x2, Val{3}())) < 10 * eps(T) / D.Δx^3
+        @test abs(derivative_right(D, x2, Val{3}())) < 10 * eps(T) / D.Δx^3
     end
 
     acc_order = 4
@@ -882,3 +876,15 @@ for T in (Float32, Float64), acc_order in (2,4,6,8)
     dest3 = D_safe*u
     @test all(i->dest1[i] ≈ dest3[i], eachindex(u))
 end
+
+
+# https://github.com/ranocha/SummationByPartsOperators.jl/issues/208
+@testset "issue #208" begin
+    for acc_order in (2, 4, 6), der_order in 1:4
+        @test_nowarn derivative_operator(Mattsson2014(),
+                                         derivative_order = der_order,
+                                         accuracy_order = acc_order,
+                                         xmin = -1.0, xmax = 1.0,
+                                         N = 20)
+    end
+end