API for Reduced Order Models using the RomApplication
The objective of this report is to present the format required to run Reduced Order Model (ROM) simulations using the RomApplication of Kratos. A simple structural example is provided to illustrate the logic of the workflow and the characteristics of the files.
The reduced order model technique currently implemented in Kratos consists of an offline-online or training-running paradigm. Therefore, these stages are presented for both, the reduced- and the hyper-reduced order models.
Make sure to install the applications required. In this example, the RomApplication and the StructuralMechanicsApplication are used. Therefore, add both these application to the Kratos configure file.
Linux:
add_app ${KRATOS_APP_DIR}/StructuralMechanicsApplication
add_app ${KRATOS_APP_DIR}/RomApplication
Windows:
CALL :add_app %KRATOS_APP_DIR%/StructuralMechanicsApplication
CALL :add_app %KRATOS_APP_DIR%/RomApplicationA dynamic structural problem using four 3D shell elements has been selected for shedding light into the expected format to use in Kratos for running a ROM simulation. The necessary files for running this example can be found here .
As can be seen in the following sketch. The model part consists of 4 shell elements and 9 nodes. Nodes 1 and 3 have constraints in their displacement degrees of freedom as well as the rotation in the z direction (perpendicular to screen). There is a point load applied in the z direction on node 9. And there are master-slave constraints between nodes 2-5 and 4-6 in the middle and 8-9 and 7-9 on the top of the part.
The expected solution for the simulation is shown next
The classical method for constructing a reduced basis is the Proper Orthogonal Decomposition (POD). This method consists on colleting a set of snapshots of the solution of a standard finite element simulation (called full order model FOM) in a matrix, to which the singular value decomposition (SVD) is applied to obtain a basis onto which the nodal finite element equations are projected.
In the MainKratos.py of the example folder, you can see that the StructuralMechanicsAnalysisMSConstraints class is a derived from StructuralMechanicsAnalysis. Inside this derived class, let us define a method to collect the nodal solution for all time steps. We do so by calling the FinalizeSolutionStep method and saving the nodal solution in a Python list. Finally, we can define a method to obtain the snapshots matrix as a Numpy array.
def FinalizeSolutionStep(self):
super().FinalizeSolutionStep()
ArrayOfDisplacements = []
for node in self._GetSolver().GetComputingModelPart().Nodes:
ArrayOfDisplacements.append(node.GetSolutionStepValue(KratosMultiphysics.ROTATION_X, 0))
ArrayOfDisplacements.append(node.GetSolutionStepValue(KratosMultiphysics.ROTATION_Y, 0))
ArrayOfDisplacements.append(node.GetSolutionStepValue(KratosMultiphysics.ROTATION_Z, 0))
ArrayOfDisplacements.append(node.GetSolutionStepValue(KratosMultiphysics.DISPLACEMENT_X, 0))
ArrayOfDisplacements.append(node.GetSolutionStepValue(KratosMultiphysics.DISPLACEMENT_Y, 0))
ArrayOfDisplacements.append(node.GetSolutionStepValue(KratosMultiphysics.DISPLACEMENT_Z, 0))
self.time_step_solution_container.append(ArrayOfDisplacements)
def GetSnapshotsMatrix(self):
SnapshotMatrix = np.zeros((len(self.time_step_solution_container[0]), len(self.time_step_solution_container)))
for i in range(len(self.time_step_solution_container)):
Snapshot_i= np.array(self.time_step_solution_container[i])
SnapshotMatrix[:,i] = Snapshot_i.transpose()
return SnapshotMatrixGiven this Numpy array containing the nodal solution of a set of simulation scenarios, the RandomizedSingularValueDecomposition returns the required basis. This basis is the truncated matrix of left singular vectors of the svd.
## Taking the SVD ###
u,s,_,_= RandomizedSingularValueDecomposition().Calculate(SnapshotMatrix,1e-6)
### Plotting singular values ###
plt.plot( range(0,len(s)), np.log(s), marker='o', markerfacecolor='blue', markersize=12, color='skyblue', linewidth=4)
plt.title('Singular Values')
plt.show()
has as many rows as degrees of freedom by the number of nodes (for this example there are 6 dofs per node: 3 displacements and 3 rotations), and it has as many columns as representative "modes" or "rom_dofs" are returned after the truncation.
Moreover, should be stored in the form of a JSON file named RomParameters.json. This JSON file must consist of "rom_settings", which in turn consists of a list with the name of the nodal degrees of freedom saved in "nodal_unknowns", and an integer stating the number of columns of the basis, stored in "number_of_rom_dofs".
Furthermore, the rows of are stored as lists of lists in "nodal_modes". The key for every list of lists is its corresponding node id.
A valid RomParameters.json for the problem at hand is shown next, in which for brevity only 2 columns of the basis have been retained. For a truncation tolerance of 1e-8, around 26 columns are kept.
{
"rom_settings": {
"nodal_unknowns": [
"ROTATION_X",
"ROTATION_Y",
"ROTATION_Z",
"DISPLACEMENT_X",
"DISPLACEMENT_Y",
"DISPLACEMENT_Z"
],
"number_of_rom_dofs": 2
},
"nodal_modes": {
"1": [
[
-0.0274537391332258,
-0.07127529747459013
],
[
-0.07875240119331663,
0.00948593116321669
],
[
1.417741611942763e-17,
1.9561040004131078e-17
],
[
6.311980879050465e-18,
7.599508201328037e-18
],
[
-1.5569829396727926e-18,
-1.0724675657707496e-17
],
[
-1.0727446748397174e-19,
-2.719138903345325e-18
]
],
"2": [
[
0.026577884352787337,
0.28393350755157803
],
[
0.120767099819073,
-0.02855243224483754
],
[
-0.0008860664136629827,
-0.0017328177630012865
],
[
0.01825515212138971,
0.0203887326444989
],
[
0.01357398473477766,
-0.026173497463154838
],
[
-0.1043712530250303,
0.2512038173812921
]
],
"3": [
[
0.02415463424401531,
0.04580523861111377
],
[
0.01849623480534586,
0.03356019538442712
],
[
0.0,
0.0
],
[
0.0,
0.0
],
[
0.0,
0.0
],
[
0.0,
0.0
]
],
"4": [
[
-0.036816510086865456,
0.05770276400767664
],
[
0.023864759349110334,
-0.07367513392524797
],
[
-0.005971855345684014,
-0.019013922634578836
],
[
0.0010817973695905184,
-0.02669476136119719
],
[
0.0007299019527745812,
-0.036813200941045114
],
[
-0.01929241731518439,
0.30451178460138395
]
],
"5": [
[
0.026577884352787327,
0.28393350755157787
],
[
0.12076709981907302,
-0.028552432244837538
],
[
-0.0008860664136629813,
-0.0017328177630012748
],
[
0.01825515212138971,
0.020388732644498843
],
[
0.01357398473477766,
-0.026173497463154845
],
[
-0.1043712530250303,
0.2512038173812921
]
],
"6": [
[
-0.036816510086865456,
0.05770276400767664
],
[
0.023864759349110334,
-0.07367513392524797
],
[
-0.005971855345684014,
-0.019013922634578836
],
[
0.0010817973695905184,
-0.02669476136119719
],
[
0.0007299019527745812,
-0.036813200941045114
],
[
-0.01929241731518439,
0.30451178460138395
]
],
"7": [
[
-0.058436773438645326,
0.08863694243164653
],
[
0.11827214091252242,
-0.04771368193062617
],
[
-0.022555472414956697,
-0.18623230608710734
],
[
0.1899753893754023,
0.3132678687137838
],
[
0.06329140162770713,
0.053066219364514854
],
[
-0.6006437762255822,
-0.3513430187068652
]
],
"8": [
[
-0.020800002079125635,
0.02693756693598514
],
[
0.0643722698434351,
-0.06363324426694524
],
[
-0.005582840070228737,
-0.12192656938106364
],
[
0.12570183262021947,
-0.185070282422351
],
[
0.02861041375552612,
-0.09130490670337461
],
[
-0.4419080755330039,
0.39509972103692226
]
],
"9": [
[
0.0,
0.0
],
[
0.0,
0.0
],
[
0.0,
0.0
],
[
0.15783861099781085,
0.06409879314571629
],
[
0.04595090769161663,
-0.019119343669429897
],
[
-0.5212759258792931,
0.021878351165028544
]
]
}
}An example on how to save in the expected format the truncated matrix of left singular vectors obtained from a structural mechanics simulation is shown next
### Saving the nodal basis ###
basis_POD={"rom_settings":{},"nodal_modes":{}}
basis_POD["rom_settings"]["nodal_unknowns"] = ["ROTATION_X","ROTATION_Y","ROTATION_Z","DISPLACEMENT_X","DISPLACEMENT_Y","DISPLACEMENT_Z"]
basis_POD["rom_settings"]["number_of_rom_dofs"] = np.shape(u)[1]
Dimensions = len(basis_POD["rom_settings"]["nodal_unknowns"])
N_nodes=np.shape(u)[0]/Dimensions
N_nodes = int(N_nodes)
node_Id=np.linspace(1,N_nodes,N_nodes)
i = 0
for j in range (0,N_nodes):
basis_POD["nodal_modes"][int(node_Id[j])] = (u[i:i+Dimensions].tolist())
i=i+Dimensions
with open('RomParameters.json', 'w') as f:
json.dump(basis_POD,f, indent=2)
print('\n\nNodal basis printed in json format\n\n')Once the RomParameters.json has been generated, one can simply call the corresponding ROM analysis. For this structural example, one should call the StructuralMechanicsAnalysisROM with the same parameters used by the original class.
Notice that even though in the training simulation the MASTER-SLAVE constraints had to be specified, this is no longer necessary because the basis already takes into account these linear constraints.
with open("ProjectParameters.json",'r') as parameter_file:
parameters = KratosMultiphysics.Parameters(parameter_file.read())
model = KratosMultiphysics.Model()
simulation = StructuralMechanicsAnalysisROM(model,parameters)
simulation.Run()The goal of the training of a hyper-reduced order model is finding a set of elements and corresponding positive weights such that by calculating the contribution to the reduced system of equations of these elements weighted by the corresponding factor, one can obtain an accurate approximation of the original reduced system, incurring a small fraction of its assembly cost.
Given a basis stored in the RomParameters.json file, one can train a Hyper Reduced Order Model by calling the appropriate solving strategy and passing to it the key "EmpiricalCubature". An example of the trainning of an HROM model for a structural mechanics simulation is given next
with open("ProjectParameters.json",'r') as parameter_file:
parameters = KratosMultiphysics.Parameters(parameter_file.read())
model = KratosMultiphysics.Model()
simulation = StructuralMechanicsAnalysisROM(model,parameters,'EmpiricalCubature')
simulation.Run()The ROM simulation will run as usual, but under the hood, the RomApplication will take care of the data management. In particular, the residual projected onto the basis will be stored for every time step. Then, the Empirical Cubature Method element selection algorithm will be called and a set of elements and weights will be saved in a json file named ElementsAndWeights.json in the folder of the case. For the example case at hand, all of 4 elements and 1 condition been selected with a weight of 1.
{
"Elements": {
"0": 1.0000000000000002,
"1": 1.0000000000000002,
"2": 1.0,
"3": 1.0
},
"Conditions": {
"0": 1.0
}
}Moreover, a hyper-reduced modelpart in MDPA format is generated. Each of the original submodelparts are now a submodelpart of the COMPUTE_HROM submodelpart. The COMPUTE_HROM submodelpart includes only the elements selected by the Empirical Cubature Method.
On the other hand, the VISUALIZE_HROM submodel part contains all of the elements on which the solution is to be printed for visualization purposes. By default, the VISUALIZE_HROM consists of all the CONDITIONS of the original modelpart. This can be costumazied by modifing the _CreateHyperReducedModelPart method in the Empirical Cubature Method class.
For running a hyper-reduced order model, it is necessary to derive a class from the appropriate ROM analysis, and import the weights to the selected elements. This can be done by calling the ModifyInitialGeometry method. An example is given next
class RunHROM(StructuralMechanicsAnalysisROM):
def ModifyInitialGeometry(self):
"""Here is the place where the HROM_WEIGHTS are assigned to the selected elements and conditions"""
super().ModifyInitialGeometry()
computing_model_part = self._solver.GetComputingModelPart()
## Adding the weights to the corresponding elements
with open('ElementsAndWeights.json') as f:
HR_data = json.load(f)
for key in HR_data["Elements"].keys():
computing_model_part.GetElement(int(key)+1).SetValue(romapp.HROM_WEIGHT, HR_data["Elements"][key])
for key in HR_data["Conditions"].keys():
computing_model_part.GetCondition(int(key)+1).SetValue(romapp.HROM_WEIGHT, HR_data["Conditions"][key])
if __name__=='__main__':
with open("ProjectParameters_HROM.json",'r') as parameter_file:
parameters = KratosMultiphysics.Parameters(parameter_file.read())
model = KratosMultiphysics.Model()
simulation = RunHROM(model,parameters)
simulation.Run()The expected format of the files required for running both ROM and HROM simulations has been presented. In summary, it can be said that a ROM simulation requires:
- A basis stored in a file named RomParameters.json.
- Calling the correct analysis strategy (for example StructuralMechanicsAnalysisROM)
On top of this, an HROM simulation also requires:
- The set of selected elements and weights, which can be stored in a file named ElementsAndWeights.json. (the weights should be imported in a derived class)
- A hyper-reduced modelpart, which only includes the necessary elements for calculation and plotting (of course, one should slightly modify the ProjectParameter.json to load this hyper-reduced model part, as well as pre-pending the COMPUTE_HROM submodel part to all the applied conditions)
