use_custom.md

-## Custom scripts
-Users can customize their use of SDPM by writing their own scripts. Each script follows a similar pattern:
-1. Creating the mesh
-2. Defining the problem
-3. Solving the flow
-4. Post-processing the results
-
-_Note_ SDPM can be used in a development or in a production environment. For the former, the code needs to be compiled, and then launched using [run.py](https://gitlab.uliege.be/am-dept/sdpm/blob/master/run.py). For the latter, the code needs to be installed to a location which is in the python path (the default install path is the user's python site-packages).
-
-### Creating the mesh
-The mesh is created using
-```python
-msh = sdpm.gmsh.MeshLoader(model_path).run(pars)
-```
-where `model_path` is the path to the file containing the model and `pars` is an optional dictionary of parameters to configure Gmsh. The model can either be a geometry (.geo) or a mesh (.msh). In the latter case, any parameters are ignored.
-
-### Defining the problem
-The problem is defined as follows
-```python
-pbl = sdpm.Problem(msh)
-pbl.setGeometry(sref, cref, xref, yref, zref, ysym)
-pbl.setFreestream(aoa, aos, minf)
-pbl.setUnsteady(n_modes, kred)
-```
-where `sref` is the reference surface area, `cref` is the reference chord length, `xref`, `yref` and `zref` are the coordinates of the reference center around which aerodynamic moments are computed, `ysym` indicates whether only a half model is provided and y-symmetry boundary condition should be applied, `aoa` and `aos` are the angles of attack and of side slip, `minf` is the freestream Mach number, `n_modes` is the number of unsteady motions and `kred` is an array of reduced frequencies. Note that the last line can be ignored for steady problems.  
-Lifting bodies are then added using
-```python
-wing = sdpm.Body(msh, body_name, te_name, cref, n_div, n_modes, n_kred)
-pbl.addBody(wing)
-```
-where `body_name` and `te_name` are the names of the lifting body to be added and of its trailing edge in the Gmsh model, `n_div` is the number of chordwise panels, and `n_kred` is the number of reduced frequencies. For steady problems, the last two arguments default to 0 and `n_div` can be set to 1.  
-Unsteady rigid body motions can be added using
-```python
-wing.addMotion()
-wing.setRigidMotion(i_mode, r_amp, d_amp, x_ref, z_ref)
-```
-where `i_mode` is the number of the motion being configured, `r_amp` is the amplitude of rotation about y, `d_amp` is the amplitude of displacement along z, and `x_ref` and `z_ref` are the x and z coordinates of the points around which the body rotates.  
-Unsteady flexible motions (mode shapes) can be added using
-```python
-x_modes, amp_modes = sdpm.utils.interpolate_modes(modes_path, wing.getNodes())
-for i in range(n_modes):
-    wing.addMotion()
-    wing.setFlexibleMotion(i, 1 / amp_modes[i], x_modes[:, 3*i:3*i+3])
-```
-where `modes_path` is the path to the file containing the mode shapes. Each row of this file must contain the x, y and z coordinates of each point on the body surface, as well as the components corresponding to the displacement along z and to the rotations around x and y. The data must then be organized as: x, y, z, dz0, rx0, ry0, dz1, rx1, ry1, etc. The utility `interpolate_modes` uses SciPy RBF implementation to interpolate the modes on the aerodynamic grid.  
-A transonic correction can be added using
-```python
-dcp = sdpm.utils.interpolate_pressure(pres_path, wing.getNodes())
-wing.setTransonicPressureGrad(dcp[:,0])
-```
-where `pres_path` is the path to the file containing the reference pressure derivative with respect to the angle of attack. Each row of this file must contain the x, y and z coordinates of each point on the body surface, as well as its associated pressure derivative. The utility `interpolate_pressure` uses SciPy RBF implementation to interpolate the pressure derivative on the aerodynamic grid.
-
-### Solving the flow
-Once the problem has been fully defined, the solver can be instantiated, updated and run using
-```python
-sol = sdpm.Solver(pbl)
-sol.update()
-sol.run()
-```
-If a transonic correction is to be applied, `SolverTransonic` can be used instead of `Solver`. The solution can be saved in Gmsh format using
-```python
-sol.save(sdpm.GmshExport(msh))
-```
-If the gradients of the solution are needed, automatic differentiation must be enabled and the solver must be run using
-```python
-adj = sdpm.Adjoint(sol)
-adj.solve()
-```
-
-### Post-processing the solution
-The solution can be accessed in Python, mainly by using
-```python
-# Steady and unsteady pressure
-sol.getPressure()
-sol.getUnsteadyPressureReal(i_mode, i_freq)
-sol.getUnsteadyPressureImag(i_mode, i_freq)
-# Steady and unsteady aerodynamic load coefficients
-sol.getLiftCoeff()
-sol.getDragCoeff()
-sol.getSideforceCoeff()
-sol.getMomentCoeff()
-sol.getUnsteadyLiftCoeffReal(i_mode, i_freq)
-sol.getUnsteadyLiftCoeffImag(i_mode, i_freq)
-sol.getUnsteadyDragCoeffReal(i_mode, i_freq)
-sol.getUnsteadyDragCoeffImag(i_mode, i_freq)
-sol.getUnsteadySideforceCoeffReal(i_mode, i_freq)
-sol.getUnsteadySideforceCoeffImag(i_mode, i_freq)
-sol.getUnsteadyMomentCoeffReal(i_mode, i_freq)
-sol.getUnsteadyMomentCoeffImag(i_mode, i_freq)
-# Generalized aerodynamic forces matrix
-sol.getGafMatrixRe(i_freq)
-sol.getGafMatrixIm(i_freq)
-```
-where `i_mode`, `i_freq` and the indices of the requested mode and frequency.  
-Additionally, chordwise pressure distributions can be extracted along the wing span by using
-```python
-sdpm.utils.create_slices(wing, vars, model_name, y_sec)
-```
-where `vars` is a dictionary containing the name and its corresponding data array (e.g. `{'cp': sol.getPressure()}`), `model_name` is the name of the model (inherited from the model file) and `y_sec` is an array of y coordinates.  
-The sensitivities of the solution can be computed using
-```python
-sens = adj.compute(of, seeds, wrt)
-```
-where `of` is the name of the variable for which the sensitivities are requested, `seeds` is an array of seed values and `wrt` is an array of names of the variables which will be differentiated with respect to. The sensitivities will be stored in the dictionary `sens` and can be accessed using `sens[wrt_var][i]` where `wrt_var` is one of the variable of `wrt` and `i` is the index of interest if the variable `wrt_var` is an array.
-
-_Note_ currently, the variables `cl`, `cd`, `cl1_re`, `cl1_im`, `cd1_re` and `cd1_im` can be differentiated with respect to `aoa` and `x_wing`.  
-_Note_ here, the sensitivity is the (adjoint) transposed matrix-vector product between the gradient and the seed. If `df` and `dx` are the seed and the sensitivity vectors, it can be expressed as  
-```python
-[dx0] = [df0_dx0 df1_dx0] [df0]  
-[dx1]   [df0_dx1 df1_dx1] [df1]
+## Custom scripts
+Users can customize their use of SDPM by writing their own scripts. Each script follows a similar pattern:
+1. Creating the mesh
+2. Defining the problem
+3. Solving the flow
+4. Post-processing the results
+
+_Note_ SDPM can be used in a development or in a production environment. For the former, the code needs to be compiled, and then launched using [run.py](https://gitlab.uliege.be/am-dept/sdpm/blob/master/run.py). For the latter, the code needs to be installed to a location which is in the python path (the default install path is the user's python site-packages).
+
+### Creating the mesh
+The mesh is created using
+```python
+msh = sdpm.gmsh.MeshLoader(model_path).run(pars)
+```
+where `model_path` is the path to the file containing the model and `pars` is an optional dictionary of parameters to configure Gmsh. The model can either be a geometry (.geo) or a mesh (.msh). In the latter case, any parameters are ignored.
+
+### Defining the problem
+The problem is defined as follows
+```python
+pbl = sdpm.Problem(msh)
+pbl.setGeometry(sref, cref, xref, yref, zref, ysym)
+pbl.setFreestream(aoa, aos, minf)
+pbl.setUnsteady(n_modes, kred)
+```
+where `sref` is the reference surface area, `cref` is the reference chord length, `xref`, `yref` and `zref` are the coordinates of the reference center around which aerodynamic moments are computed, `ysym` indicates whether only a half model is provided and y-symmetry boundary condition should be applied, `aoa` and `aos` are the angles of attack and of side slip, `minf` is the freestream Mach number, `n_modes` is the number of unsteady motions and `kred` is an array of reduced frequencies. Note that the last line can be ignored for steady problems.  
+Lifting bodies are then added using
+```python
+wing = sdpm.Body(msh, body_name, te_name, cref, n_div, n_modes, n_kred)
+pbl.addBody(wing)
+```
+where `body_name` and `te_name` are the names of the lifting body to be added and of its trailing edge in the Gmsh model, `n_div` is the number of chordwise panels, and `n_kred` is the number of reduced frequencies. For steady problems, the last two arguments default to 0 and `n_div` can be set to 1.  
+Unsteady rigid body motions can be added using
+```python
+wing.addMotion()
+wing.setRigidMotion(i_mode, r_amp, d_amp, x_ref, z_ref)
+```
+where `i_mode` is the number of the motion being configured, `r_amp` is the amplitude of rotation about y, `d_amp` is the amplitude of displacement along z, and `x_ref` and `z_ref` are the x and z coordinates of the points around which the body rotates.  
+Unsteady flexible motions (mode shapes) can be added using
+```python
+x_modes, amp_modes = sdpm.utils.interpolate_modes(modes_path, wing.getNodes())
+for i in range(n_modes):
+    wing.addMotion()
+    wing.setFlexibleMotion(i, 1 / amp_modes[i], x_modes[:, 3*i:3*i+3])
+```
+where `modes_path` is the path to the file containing the mode shapes. Each row of this file must contain the x, y and z coordinates of each point on the body surface, as well as the components corresponding to the displacement along z and to the rotations around x and y. The data must then be organized as: x, y, z, dz0, rx0, ry0, dz1, rx1, ry1, etc. The utility `interpolate_modes` uses SciPy RBF implementation to interpolate the modes on the aerodynamic grid.  
+A transonic correction can be added using
+```python
+dcp = sdpm.utils.interpolate_pressure(pres_path, wing.getNodes())
+wing.setTransonicPressureGrad(dcp[:,0])
+```
+where `pres_path` is the path to the file containing the reference pressure derivative with respect to the angle of attack. Each row of this file must contain the x, y and z coordinates of each point on the body surface, as well as its associated pressure derivative. The utility `interpolate_pressure` uses SciPy RBF implementation to interpolate the pressure derivative on the aerodynamic grid.
+
+### Solving the flow
+Once the problem has been fully defined, the solver can be instantiated, updated and run using
+```python
+sol = sdpm.Solver(pbl)
+sol.update()
+sol.run()
+```
+If a transonic correction is to be applied, `SolverTransonic` can be used instead of `Solver`. The solution can be saved in Gmsh format using
+```python
+sol.save(sdpm.GmshExport(msh))
+```
+If the gradients of the solution are needed, automatic differentiation must be enabled and the solver must be run using
+```python
+adj = sdpm.Adjoint(sol)
+adj.solve()
+```
+
+### Post-processing the solution
+The solution can be accessed in Python, mainly by using
+```python
+# Steady and unsteady pressure
+sol.getPressure()
+sol.getUnsteadyPressureReal(i_mode, i_freq)
+sol.getUnsteadyPressureImag(i_mode, i_freq)
+# Steady and unsteady aerodynamic load coefficients
+sol.getLiftCoeff()
+sol.getDragCoeff()
+sol.getSideforceCoeff()
+sol.getMomentCoeff()
+sol.getUnsteadyLiftCoeffReal(i_mode, i_freq)
+sol.getUnsteadyLiftCoeffImag(i_mode, i_freq)
+sol.getUnsteadyDragCoeffReal(i_mode, i_freq)
+sol.getUnsteadyDragCoeffImag(i_mode, i_freq)
+sol.getUnsteadySideforceCoeffReal(i_mode, i_freq)
+sol.getUnsteadySideforceCoeffImag(i_mode, i_freq)
+sol.getUnsteadyMomentCoeffReal(i_mode, i_freq)
+sol.getUnsteadyMomentCoeffImag(i_mode, i_freq)
+# Generalized aerodynamic forces matrix
+sol.getGafMatrixRe(i_freq)
+sol.getGafMatrixIm(i_freq)
+```
+where `i_mode`, `i_freq` and the indices of the requested mode and frequency.  
+Additionally, chordwise pressure distributions can be extracted along the wing span by using
+```python
+sdpm.utils.create_slices(wing, vars, model_name, y_sec)
+```
+where `vars` is a dictionary containing the name and its corresponding data array (e.g. `{'cp': sol.getPressure()}`), `model_name` is the name of the model (inherited from the model file) and `y_sec` is an array of y coordinates.  
+The derivatives of the Generalized Aerodynamic Force (GAF) matrices can be calculated and accessed using
+```python
+grd = sdpm.Gradient(sol)
+grd.run()
+grd.computePartialsGafRe(i_freq, i_mode, i_row, j_col, wrt)
+grd.computePartialsGafIm(i_freq, i_mode, i_row, j_col, wrt)
+```
+The last two lines will return the derivative of the real and imaginary parts of the GAF matrix entry at row `i_row` and column `j_col` with respect to the component `wrt` of mode number `i_mode` at frequency `i_freq` for all body nodes. `wrt` can either be `'a'` (rigid pitch), `'h'` (rigid plunge), `'dz'` (modal displacements along z), `'rx'` (modal rotations about x) or `'ry'` (modal rotations about y).  
+If automatic differentiation has been enabled, the sensitivities of the solution can be computed using
+```python
+sens = adj.compute(of, seeds, wrt)
+```
+where `of` is the name of the variable for which the sensitivities are requested, `seeds` is an array of seed values and `wrt` is an array of names of the variables which will be differentiated with respect to. The sensitivities will be stored in the dictionary `sens` and can be accessed using `sens[wrt_var][i]` where `wrt_var` is one of the variable of `wrt` and `i` is the index of interest if the variable `wrt_var` is an array.
+
+_Note_ currently, the variables `cl`, `cd`, `cl1_re`, `cl1_im`, `cd1_re` and `cd1_im` can be differentiated with respect to `aoa` and `x_wing`.  
+_Note_ here, the sensitivity is the (adjoint) transposed matrix-vector product between the gradient and the seed. If `df` and `dx` are the seed and the sensitivity vectors, it can be expressed as  
+```python
+[dx0] = [df0_dx0 df1_dx0] [df0]  
+[dx1]   [df0_dx1 df1_dx1] [df1]
 ```
\ No newline at end of file