diff --git a/src/SpinShuttling.jl b/src/SpinShuttling.jl
index 380bdee..cc5ab7a 100644
--- a/src/SpinShuttling.jl
+++ b/src/SpinShuttling.jl
@@ -197,7 +197,7 @@ end
 """
 Calculate the dephasing matrix of a given spin shuttling model.
 """
-function dephasingmatrix(model::ShuttlingModel)::Symmetric{<:Real}
+function dephasingmatrix(model::ShuttlingModel)::Matrix{<:Real}
     n = model.n
     W = zeros(2^n, 2^n)
     for j in 1:2^n
@@ -208,7 +208,7 @@ function dephasingmatrix(model::ShuttlingModel)::Symmetric{<:Real}
             W[k, j] = W[j, k]
         end
     end
-    return Symmetric(W)
+    return W
 end
 
 """
diff --git a/test/runtests.jl b/test/runtests.jl
index f761b51..909b97c 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -4,8 +4,8 @@ using SpinShuttling
 using UnicodePlots
 using ProgressMeter
 
-include("testintegration.jl")
-include("testfidelity.jl")
-include("teststochastics.jl")
-include("testspectrum.jl")
+# include("testintegration.jl")
+# include("testfidelity.jl")
+# include("teststochastics.jl")
+# include("testspectrum.jl")
 include("testdephasingfactor.jl")
\ No newline at end of file
diff --git a/test/testdephasingfactor.jl b/test/testdephasingfactor.jl
index c8688c2..1158cd2 100644
--- a/test/testdephasingfactor.jl
+++ b/test/testdephasingfactor.jl
@@ -1,46 +1,59 @@
 ## 
-@testset begin "test dephasing matrix"
-    T=400; L=10; σ = sqrt(2) / 20; M = 20000; N=601; κₜ=1/20;κₓ=1/0.1;
-
-    B=OrnsteinUhlenbeckField(0,[κₜ,κₓ],σ)
-    model=OneSpinModel(T,L,N,B)
+@testset begin
+    "test dephasing matrix"
+    T = 400
+    L = 10
+    σ = sqrt(2) / 20
+    M = 20000
+    N = 601
+    κₜ = 1 / 20
+    κₓ = 1 / 0.1
+
+    B = OrnsteinUhlenbeckField(0, [κₜ, κₓ], σ)
+    model = OneSpinModel(T, L, N, B)
     # test customize println
     println(model)
 
-    f=statefidelity(model)
-    w=dephasingmatrix(model)
-
-    w2=SpinShuttling.parallelsampling(i->dephasingmatrix(model,randn(N)),M)[1]
-    @test norm(w-w2)<2e-2
-    w3=sampling(i->dephasingmatrix(model,randn(N)),M)[1]
-    @test norm(w-w3)<2e-2
+    f = statefidelity(model)
+    w = dephasingmatrix(model)
 
+    w2 = sampling(model, dephasingmatrix, M)
+    @test norm(w - w2) < 2e-2
 
-    rho=model.Ψ*model.Ψ'
-    @test rho[1,1]+rho[2,2]≈ 1
-    @test w==w'
+    rho = model.Ψ * model.Ψ'
+    @test rho[1, 1] + rho[2, 2] ≈ 1
+    @test w == w'
 
     println(w)
-    f_c=(model.Ψ'*(w.*rho)*model.Ψ)
-
-    @test f≈ f_c
+    f_c = (model.Ψ' * (w .* rho) * model.Ψ)
+    f_s = (model.Ψ' * (abs.(w2) .* rho) * model.Ψ)
+    @test f ≈ f_c
+    @test f ≈ f_s
 end
 
 ## 
-@testset begin "test two spin dephasing matrix"
-    L=10; σ =sqrt(2)/20; M=20000; N=501; T1=200; T0=25*0.05; κₜ=1/20; κₓ=1/0.1;
-    B=OrnsteinUhlenbeckField(0,[κₜ,κₓ],σ)
-    model=TwoSpinModel(T0, T1, L, N, B)
+@testset begin
+    "test two spin dephasing matrix"
+    L = 10
+    σ = sqrt(2) / 20
+    M = 20000
+    N = 501
+    T1 = 200
+    T0 = 25 * 0.05
+    κₜ = 1 / 20
+    κₓ = 1 / 0.1
+    B = OrnsteinUhlenbeckField(0, [κₜ, κₓ], σ)
+    model = TwoSpinModel(T0, T1, L, N, B)
     println(model)
-    f=statefidelity(model)
-    w=dephasingmatrix(model)
+    f = statefidelity(model)
+    w = dephasingmatrix(model)
 
-    rho=model.Ψ*model.Ψ'
-    @test sum([rho[i,i] for i in 1:4 ]) ≈ 1
-    @test w==w'
+    rho = model.Ψ * model.Ψ'
+    @test sum([rho[i, i] for i in 1:4]) ≈ 1
+    @test w == w'
 
     println(w)
-    f_c=(model.Ψ'*(w.*rho)*model.Ψ)
+    f_c = (model.Ψ' * (w .* rho) * model.Ψ)
 
     println(f)
     println(f_c)