diff --git a/Dipolar-Gas-Simulator/+VariationalSolver2D/@Calculator/calculateVariationalEnergy.m b/Dipolar-Gas-Simulator/+VariationalSolver2D/@Calculator/calculateVariationalEnergy.m
index f3466ae..cd16eec 100644
--- a/Dipolar-Gas-Simulator/+VariationalSolver2D/@Calculator/calculateVariationalEnergy.m
+++ b/Dipolar-Gas-Simulator/+VariationalSolver2D/@Calculator/calculateVariationalEnergy.m
@@ -1,4 +1,4 @@
-function E = calculateVariationalEnergy(~, psi, Params, ell, Transf, VDk, V)
+function E = calculateVariationalEnergy(~, psi, Params, ell, Transf, V)
                                
     g_eff       = Params.gs      * (1/(sqrt(2*pi)*ell));
     gamma_eff   = Params.gammaQF * (sqrt(2/5)/(pi^(3/4)*ell^(3/2)));
@@ -9,6 +9,8 @@ function E = calculateVariationalEnergy(~, psi, Params, ell, Transf, VDk, V)
     normfac     = Params.Lx*Params.Ly/numel(psi);
     
     % DDIs
+    VParams.ell = ell;
+    VDk         = this.Calculator.calculateVDkWithCutoff(Transf, Params, VParams); % VDk depends on the variational parameters. THIS HAS TO BE RECALCULATED HERE FOR THE CONSTRAINED OPTIMIZATION TO WORK!
     frho        = fftn(abs(psi).^2);
     Phi         = real(ifftn(frho.*VDk));
     Eddi        = 0.5*Params.gdd*Phi.*abs(psi).^2/(sqrt(2*pi)*ell);%
diff --git a/Dipolar-Gas-Simulator/+VariationalSolver2D/@DipolarGas/initialize.m b/Dipolar-Gas-Simulator/+VariationalSolver2D/@DipolarGas/initialize.m
index b93ccae..81d22b2 100644
--- a/Dipolar-Gas-Simulator/+VariationalSolver2D/@DipolarGas/initialize.m
+++ b/Dipolar-Gas-Simulator/+VariationalSolver2D/@DipolarGas/initialize.m
@@ -5,8 +5,6 @@ function [psi,V,VDk] = initialize(this,Params,VParams,Transf)
     
     % == Calculating the DDIs == %
     VDk = this.Calculator.calculateVDkWithCutoff(Transf, Params, VParams);
-    save(sprintf(strcat(this.SaveDirectory, '/VDk_M.mat')),'VDk');
-    fprintf('Computed and saved DDI potential in Fourier space with cutoff.\n')
     
     % == Setting up the initial wavefunction == %
     psi = this.setupWavefunction(Params,Transf);
diff --git a/Dipolar-Gas-Simulator/+VariationalSolver2D/@DipolarGas/run.m b/Dipolar-Gas-Simulator/+VariationalSolver2D/@DipolarGas/run.m
index e250f0c..987e2f1 100644
--- a/Dipolar-Gas-Simulator/+VariationalSolver2D/@DipolarGas/run.m
+++ b/Dipolar-Gas-Simulator/+VariationalSolver2D/@DipolarGas/run.m
@@ -34,7 +34,7 @@ function [Params, Transf, psi, V, VDk] = run(this)
         
         t_idx                 = 1; % Start at t = 0;
 
-        [psi,V,VDk]           = this.initialize(Params,VParams,Transf); % Recalculate Psi, VDk with new value for the variational parameter
+        [~,V,VDk]             = this.initialize(Params,VParams,Transf); % Recalculate VDk with new value for the variational parameter, keep new psi at the end of loop to increase likelihood of faster convergence
                                                                         
         % --- Adding some noise --- 
         % Noise added in every iteration to ensure it is not stuck in some local minimum