diff --git a/docs/Project.toml b/docs/Project.toml
index b3c1ed2088..b44ae23671 100644
--- a/docs/Project.toml
+++ b/docs/Project.toml
@@ -6,6 +6,7 @@ DomainSets = "5b8099bc-c8ec-5219-889f-1d9e522a28bf"
 Flux = "587475ba-b771-5e3f-ad9e-33799f191a9c"
 Integrals = "de52edbc-65ea-441a-8357-d3a637375a31"
 IntegralsCubature = "c31f79ba-6e32-46d4-a52f-182a8ac42a54"
+LineSearches = "d3d80556-e9d4-5f37-9878-2ab0fcc64255"
 Lux = "b2108857-7c20-44ae-9111-449ecde12c47"
 ModelingToolkit = "961ee093-0014-501f-94e3-6117800e7a78"
 NeuralPDE = "315f7962-48a3-4962-8226-d0f33b1235f0"
@@ -35,6 +36,6 @@ OptimizationOptimisers = "0.1"
 OptimizationPolyalgorithms = "0.1"
 OrdinaryDiffEq = "6.31"
 Plots = "1.36"
-QuasiMonteCarlo = "0.2"
+QuasiMonteCarlo = "0.3"
 Roots = "2.0"
 SpecialFunctions = "2.1"
diff --git a/docs/make.jl b/docs/make.jl
index 67ee53b3d4..09bc53388c 100644
--- a/docs/make.jl
+++ b/docs/make.jl
@@ -4,6 +4,7 @@ cp("./docs/Manifest.toml", "./docs/src/assets/Manifest.toml", force = true)
 cp("./docs/Project.toml", "./docs/src/assets/Project.toml", force = true)
 
 ENV["GKSwstype"] = "100"
+ENV["JULIA_DEBUG"] = "Documenter"
 using Plots
 
 include("pages.jl")
diff --git a/docs/src/examples/ks.md b/docs/src/examples/ks.md
index bb2b1064f5..ae960b5c5b 100644
--- a/docs/src/examples/ks.md
+++ b/docs/src/examples/ks.md
@@ -27,7 +27,7 @@ where $\theta = t - x/2$ and with initial and boundary conditions:
 We use physics-informed neural networks.
 
 ```@example ks
-using NeuralPDE, Lux, ModelingToolkit, Optimization, OptimizationOptimJL
+using NeuralPDE, Lux, ModelingToolkit, Optimization, OptimizationOptimJL, LineSearches
 import ModelingToolkit: Interval, infimum, supremum
 
 @parameters x, t
@@ -71,7 +71,7 @@ callback = function (p, l)
     return false
 end
 
-opt = OptimizationOptimJL.BFGS()
+opt = OptimizationOptimJL.LBFGS(linesearch = BackTracking())
 res = Optimization.solve(prob, opt; callback = callback, maxiters = 2000)
 phi = discretization.phi
 ```
diff --git a/docs/src/examples/linear_parabolic.md b/docs/src/examples/linear_parabolic.md
index 0ae1432763..fa5bb8b6d7 100644
--- a/docs/src/examples/linear_parabolic.md
+++ b/docs/src/examples/linear_parabolic.md
@@ -27,6 +27,7 @@ with a physics-informed neural network.
 using NeuralPDE, Lux, ModelingToolkit, Optimization, OptimizationOptimJL
 using Plots
 import ModelingToolkit: Interval, infimum, supremum
+using LineSearches
 
 @parameters t, x
 @variables u(..), w(..)
@@ -92,7 +93,7 @@ callback = function (p, l)
     return false
 end
 
-res = Optimization.solve(prob, BFGS(); callback = callback, maxiters = 5000)
+res = Optimization.solve(prob, LBFGS(linesearch = BackTracking()); callback = callback, maxiters = 5000)
 
 phi = discretization.phi