Numerical modeling & simulation are today’s first priority for all industries which would like to set time- and cost-efficient design process. But needs are more and more complex : non linear physic simulation (turbulence, vibrations, damage propagation), fast dynamic simulation (blade loss, missile impact) and simulation on various scale (boundary layer of a wind turbine and the complete wind farm in a same time). As computation times is expensive, it is essential to find new tools in order to generate lighter surrogate model. It is the role of the Proper Orthogonal Decomposition (POD) tool that is presented below.
What is possible with this new tool
Thanks to simple methods already implemented in Scilab (Singular Value Decomposition on space/time correlation matrices), this tool generates an orthonormal basis of modes based on experimental or simulation data. Taking the first more energetics modes, depending on the case, can sum up to 99% of the total energy of the system while drastically reducing the problem size . Applying a mathematical-based model is not anymore about solving equations on each nodes of the mesh, but solving equations on each mode of the system behavior. Here lies the lighter model.
New way of simulating
Summarizing a system behavior in a few modes allows many new applications. It is now possible to make up to 10.000 times faster simulations, useful for time-critical applications such as Model Predictive Control . Thanks to multi-parametric prediction, it is also possible to make fast design space exploration, for design optimization process for instance . The perfect complementary tool for all prototyping tool.
This paper will show you the good use of the toolbox around the case of a Danielson reactor model (basically a CVD reactor using RF coil to inductively heat the substrate). This CFD study has been run at 69 sparse grid points for 3 varying parameters : Frequency (20 – 70kHz), RF current in coils (150 – 300 A) and Flow (0 – 29 slm). Let see how behaves the imaginary part of the static potential in the cavity and get more informations on non simulated parameters case. It is a steady case.
To set up the learning phase, you will first need to set up your mesh in Scilab. To do so you will need two text files that you can transform into array in Scilab using the fscanfMat function. The first one gathers X, Y nodes coordinates, one column for each. The second one is the cell description (connectivity) of your mesh. This is an array where the nth line is composed of the position of each nodes (regarding the first coordinates file) composing the nth cell.
In the example on the left, we can see the the X and Y coordinates on the left on two columns. Beside lies the connectivity file, where we can see that the mesh is composed of quad cells (4 columns) and that the first cell is composed of the 6219, 6218, 6294, 6295 nodes in the coordinates list.
Now import your field data, in our case, the imaginary part of the static potential (Az_i) field. To do so, gather the nodal values for all time/parameter simulations (one snapshot per column) and follow the same process as above.
Create mesh and simulation
It’s time to see if all the import phase (mesh + data) went right. To do so, you will need to create a Mesh and a Simulation. Those are the two main objects of the Scilab POD toolbox.
Mesh object properties : node_coords / nb_nodes / cell_type / elements / nb_elements / border_edges / border_normal_vectors / element_x / element_y / element_area
This function will allow you to create a mesh object, depending on the nodes coordinates (node_coords) and cell description (elements) you previously imported. The boolean input argument compute_data just gives your mesh more properties.
Simulation object properties : _mesh / nb_samples/ values
This function will allow you to create a simulation object, depending on the mesh object you just created and the nodal values (node_values) of your field.
Try it and validate the result with the show_simulation function. You just have to insert your simulation object (simulation) you just created and it will plot your field as on above image. Depending on the number of snapshot you have, you can display a single frame or a complete animation. If you set varargin to %t, it will save you all the frame as .gif images.
What can be useful when testing this toolbox is to save your mesh and/or simulation object in order to avoid all the previous process that can be annoying. You would then be able to load it in Scilab easily.
Generate reduced modal basis
What you get from now : A mesh, a learning simulation.
The idea behind POD is to sum up your system behavior into a few of the most energetic spatial (classical method) or temporal (snapshots method) correlations called modes . Theoretically, those modes constitute an orthonormal basis of your simulation domain. Those modes are sorted in order for you to consider only the more energetic one, and thus reduce your problem complexity. Let see how to generate this modal basis in Scilab.
Limit_conditions (generated by the prepare_limit_conditions function) is also an object in this toolbox gathering this only two properties : homogeneous (boolean) and time_invariant (boolean). This is mandatory to set a mathematical based model reduction. In this case, such accuracy is not needed. So let just say that your simulation is not time invariant (%f) and that your boundary conditions are homogeneous (%t). Doing so you will just avoid analyzing the mean flow, and get a better display of your modes.
Pod_basis object properties : u_m (mean field) / u_s (varying field) / phi (modes) / lambda (associated eigenvalues)
Here is the function that generates your modal basis … not reduced, yet. The basis generation could take a bit time because of the construction/diagonalization of a huge correlation matrix. For now, a 100.000 nodes mesh, 1000 frames/snapshots takes 12min computing on Scilab 6 to get the pod_basis.
Reduced_basis object properties : n_pod (number of selected modes) / phi (n_pod first modes) / u_m
Now that you get your pod_basis, you can reduce it keeping only the most energetic modes and visualize them (modes and associated energy). To do so, you just have to specify your simulation (simulation_plan), the pod_basis you just generated (pod_basis) and the number of modes you want to keep (varargin, default 12). You will get something like the above image.
Be careful ! As said earlier, the pod basis is theoretically an orthonormal basis, meaning that keeping the full basis will give you the exact solution you would have simulated. But, due to some numerical error, taking too much modes into account increases the reconstruction/prediction error. Take the time to look for the optimal number of modes regarding accuracy vs computing time.
I you want to evaluate the accuracy of your reduced basis, or may be compress your data, try to project your simulation on this new basis.
Reconstruct_simulation is a simple simulation object based on a simple one field frame (node_values) projection on the reduced basis (reduced_basis). Besides making the projection, the reconstruct_simulation function will display the n_pod projection coefficients. Just use the show_simulation function in order to see the simulation reconstruction.
The compute_rel_error function returns the normalized root mean square deviation between your nodal value of field (column vector expected) and the reconstruct_simulation nodes values (value) you get at the previous stage.
In the above image, you can see the reconstruction of the Az_i field for the 69th parametric case and, next to it, the relative error computed for all the parametric cases. For 12 modes, the maximum error reach 4%.
Making time/parametric prediction
What you get from now : A mesh, a learning simulation, a reduced pod basis.
The idea is that, with a few simulations on specific parameters, you can sum up the behavior of your system in any conditions. This can take a Design of Experiment to set up relevant simulation cases. But once this is done, and the associated reduced basis generated following the previous process, you can now predict your system behavior on non-simulated cases.
Practically it consists in interpolating the projection coefficients following specific methods. But the appropriate interpolation method depends on your simulation, because your coefficients are time- and parameter-dependent. Here are the two methods currently available in the toolbox.
Steady multi-parametric interpolation
This function will first compute the projection coefficients on the reduced basis (basis – reduced_basis object) for all of your learning simulations (lrnPlan – simulation object) and will interpolate the coefficients (NewPInterp – array) for the new parameter (NewParam – array) with Radial Basis Function . To do so, you just need to set your entire parameter space grid (DOE – array) composed of your parameters values for each simulation case, one parameter per column.
Unsteady mono-parametric interpolation
Unsteady simulation are more complex to predict because we can’t make a simple interpolation on the coefficients, considering that those are time- and parameter-dependent. So you will have to make a second POD on the reduced basis you got, using n2_pod function .
Giving the number of learning parameters (np), the first reduced_basis (reduced_bas) and the learning simulation (training_simu), this function returns a new POD basis (podi) and the associated coefficients (cofi). Those coefficients are used for the interpolation in the following function :
Just as RBF_interpolation, this function will return the interpolated coefficients (NewPInterp) for the new value of parameter (NewParam), given the parameter space grid (DOE – have to be line vector), the number of modes n_pod in the first reduced basis (NFstMod) and the coefficients given by the n2_pod function (cof).
Making prediction from NewPInterp
What you get from now : A mesh, learning simulation, reduced pod basis and interpolated coefficients NewPInterp.
The predict_simulation function creates the simulation (predicted_simulation) you want to see. Input arguments are simply the NewPInterp coefficients you just got (Coef), the [initial] reduced_basis (reduced_basis), and the mesh (_mesh). If you have used the n2_pod function to get interpolated coefficients for your unsteady simulation, set n2 input argument to %t and add the second POD basis podi generated by n2_pod function (varargin).
Want to see your prediction ? Just use show_simulation function!
ODE based model order reduction
 Huisman, L. (2005). Control of glass melting processes based on reduced CFD models, Eindhoven : Technische Universiteit Eindhoven D OI : 10.6100/IR586128
 Hovland S., Willcox K. and Gravdahl J.T., MPC for large-scale sone or no input.ystems via model reduction and multiparametric quadratic programming.
 Patrick Allen LeGresley, Application of proper orthogonal decomposition to design decomposition methods, October 2005.
 Walton S., Hassan O., Morgan K., Reduced order modelling for unsteady fluid flow using proper orthogonal decomposition and radial basis functions, College of Engineering, Swansea University, 2013.