diff --git a/Project.toml b/Project.toml
index dbc76e2..72dc901 100644
--- a/Project.toml
+++ b/Project.toml
@@ -19,7 +19,6 @@ SciMLBase = "0bca4576-84f4-4d90-8ffe-ffa030f20462"
 SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
 StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91"
 Triangulate = "f7e6ffb2-c36d-4f8f-a77e-16e897189344"
-UnPack = "3a884ed6-31ef-47d7-9d2a-63182c4928ed"
 
 [extras]
 Enzyme = "7da242da-08ed-463a-9acd-ee780be4f1d9"
diff --git a/test/testoptimization.jl b/test/testoptimization.jl
index 0a8eb9a..3b8a4e9 100755
--- a/test/testoptimization.jl
+++ b/test/testoptimization.jl
@@ -1,12 +1,12 @@
-module enzymeDiff
 
 using dJUICE
 using MAT
 using Test
 using Enzyme
-
+Enzyme.API.typeWarning!(false)
 Enzyme.Compiler.RunAttributor[] = false
-using Optimization, OptimizationOptimJL
+
+using Optimization
 
 #Load model from MATLAB file
 #file = matopen(joinpath(@__DIR__, "..", "data","temp12k.mat")) #BIG model
@@ -26,21 +26,32 @@ md.inversion.max_parameters = ones(md.mesh.numberofvertices)*(1.0e3)
 md.inversion.independent_string = "FrictionCoefficient"
 
 α = md.inversion.independent
+α[2] = 10.0
 ∂J_∂α = zero(α)
 
 femmodel=dJUICE.ModelProcessor(md, :StressbalanceSolution)
 n = length(α)
 # use user defined grad, errors!
-optprob = OptimizationFunction(dJUICE.costfunction, Optimization.AutoEnzyme(), grad=dJUICE.computeGradient)
-#optprob = OptimizationFunction(dJUICE.costfunction, Optimization.AutoEnzyme())
-#prob = Optimization.OptimizationProblem(optprob, α, femmodel, lb=md.inversion.min_parameters, ub=md.inversion.max_parameters)
-prob = Optimization.OptimizationProblem(optprob, α, femmodel)
-sol = Optimization.solve(prob, Optim.LBFGS())
+optprob = OptimizationFunction(DJUICE.costfunction, Optimization.AutoEnzyme())
+prob = Optimization.OptimizationProblem(optprob, α, femmodel, lb=md.inversion.min_parameters, ub=md.inversion.max_parameters, maxiters=1000)
+#prob = Optimization.OptimizationProblem(optprob, α, femmodel)
+#sol = Optimization.solve(prob, Optim.NelderMead())
+sol = Optimization.solve(prob, Optimization.LBFGS())
+
+# compute the gradient by enzyme
+md = DJUICE.solve(md, :grad)
+enz_grad = md.results["StressbalanceSolution"]["Gradient"]
+
+# evaluating the gradient
+sol.cache.f.grad(∂J_∂α, prob.u0, prob.p)
 
+#independent_enum = DJUICE.StringToEnum(md.inversion.independent_string)
+#DJUICE.InputUpdateFromVectorx(femmodel, sol.u, independent_enum, DJUICE.VertexSIdEnum)
+#DJUICE.RequestedOutputsx(femmodel, [independent_enum])
+#DJUICE.OutputResultsx(femmodel, md, :StressbalanceSolution)
 
 #md = solve(md, :sb)
 
 # compute gradient by finite differences at each node
 #addJ = md.results["StressbalanceSolution"]["Gradient"]
 
-end