Constraining-densities-and-potentials.md

-In each iteration used can impose *by hand* various constraints on densities or potentials. [problem-definition.h](https://gitlab.fizyka.pw.edu.pl/wtools/wslda/-/blob/public/st-project-template/problem-definition.h) file provides handlers for these operations.  
+[[_TOC_]]
+# Imposing constraints on densities and potentials
+In each iteration, the user can impose *by hand* various constraints on densities or potentials. [problem-definition.h](https://gitlab.fizyka.pw.edu.pl/wtools/wslda/-/blob/public/st-project-template/problem-definition.h) file provides handlers for these operations.  
 
-If you want to constrain density you need to provide a body of the function:
+If you want to constrain density, you need to provide a body of the function:
 ```c
 /**
  * THIS FUNCTION IS CALLED DURING THE SELF-CONSISTENT PROCESS.
@@ -36,11 +38,11 @@ void modify_densities(int it, wslda_density h_densities, double *params, size_t
 ```
   
 
-If you want to constrain potentials you need to provide a body of the function:
+If you want to constrain potentials, you need to provide a body of the function:
 ```c
 /**
  * THIS FUNCTION IS CALLED DURING THE SELF-CONSISTENT PROCESS.
- * Before each diagonalization process, user can modify arbitrarily potentials
+ * Before each diagonalization process, user can arbitrarily modify potentials
  * @param it iteration number
  * @param h_densities structure with densities, see (wiki) documentation for list of fields
  *                    NOTE: densities structure is processed by modify_densities(...) function before call ot this function.
@@ -74,8 +76,8 @@ void modify_potentials(int it, wslda_density h_densities, wslda_potential h_pote
 
 ```
 
-# Ordering of function calls
-Function are called as follow:
+## Order of function calls
+Functions are called as follows:
 ```c
 while( not (convergence criteria) )
 {
@@ -90,12 +92,12 @@ while( not (convergence criteria) )
 }
 ```
 
-# Example: imprinting vortex in the center of the system
-The standard method of imprinting quantum vortex is to impose that order parameter has the following structure:  
+## Example: imprinting vortex in the center of the system
+The standard method of imprinting quantum vortex is to impose that the order parameter has the following structure:  
 ```math
 \Delta(\vec{r})=|\Delta(\vec{r})|e^{i\phi}
 ```
-where $`\phi`$ is angle in cylindrical system, i.e $`\phi=\arctan(y/x)`$. The procedure is implemented in function below:
+where $`\phi`$ is angle in cylindrical system, i.e $`\phi=\arctan(y/x)`$. The procedure is implemented in the function below:
 
 
 ```c
@@ -121,4 +123,22 @@ void modify_potentials(int it, wslda_density h_densities, wslda_potential h_pote
         ixyz++; // go to next point,  it should be last line of the triple loop
     }
 }
-```
+```
\ No newline at end of file
+
+# Constrained density functional theory calculations
+In the standard approach, we provide an external potential $`V_{ext}(\bm{r})`$, and after the minimization process, we obtain the density distribution. However, we can treat the external potential $`V_{ext}(\bm{r})`$ as a generalized Lagrange multiplier that controls the density. For example, we can ask: _What is the shape of the external potential that provides the desired density distribution_ $`n_0(\bm{r})`$? Such a type of minimization is popular in nuclear physics and can also be implemented with API provided by W-SLDA Toolkit.
+In order to find such a potential, one can introduce a procedure similar to adjusting the chemical potential:
+```math
+V_{ext}^{(it+1)}(\bm{r}) := V_{ext}^{(it)}(\bm{r}) + ab^{it}\frac{n^{(it)}(\bm{r}) - n_0(\bm{r})}{n_0(\bm{r}) + c}
+``` 
+where $`a,b,c`$ are parameters that need to be adjusted to achieve convergence, $`it`$ stands for the iteration number. 
+
+As an example, let us consider the quasi-1D case, where search for the external potential that produces density distribution in the form of two Gaussians:
+```math
+n_0(\bm{r}) = A_1\exp\left(-\frac{(x-x_1)^2}{2\sigma_1^2}\right) + A_2\exp\left(-\frac{(x-x_2)^2}{2\sigma_2^2}\right)
+```
+In the figure below, the requested density is plotted by red, while with the blue color, we plot the density obtained after iterating the code with the prescription for the update of the external potential $`V_{ext}^{(it+1)}(\bm{r})`$:
+![image](uploads/89ebc8d36f87f9f2000f0a17e0701c73/image.png)
+
+The resulting external potential has a form as presented in the figure below:
+![image](uploads/dd648890418d2b9c87db7a2c7a697c17/image.png)
\ No newline at end of file