Généralités
Les études paramétriques sont un moyen utile et efficace pour adapter des simulations. Les objectifs peuvent être de nature différente : Trouver une structure optimale avec un poids minimal tout en garantissant la sécurité et en économisant des ressources, adapter la finesse du maillage pour obtenir des résultats de simulation fiables ou faire varier les propriétés de matériau pour vérifier la sensibilité de la simulation et pouvoir faire des déclarations fiables sur la sécurité structurelle. Pour ces cas et bien d'autres, une modélisation paramétrée est appropriée. Pour ne pas avoir à calculer toutes les variantes individuellement, l'API Dlubal peut être utilisée. Cela permet d'exécuter le calcul de l'étude paramétrique sans avoir à le superviser.
Modélisation paramétrée
De nombreuses fonctions dans RFEM 6, RSTAB 9 et RSECTION sont déjà paramétrables. De plus amples informations sont disponibles dans les manuels correspondants :
- Manuels en ligne RFEM 6 | Fonctions du programme | Saisie paramétrée
- Manuels en ligne RSTAB 9 | Fonctionnalité de produit | Saisie paramétrée
- Manuels en ligne de RSECTION 1 | Fonctionnalités du produit | Saisie paramétrée
Accès au modèle via l'API
Pour pouvoir accéder au modèle via l'API, le package correspondant doit être installé. Vous pouvez également utiliser directement la console ou le gestionnaire de scripts dans le programme principal au lieu d'un éditeur externe. De plus amples informations sont disponibles sous les liens suivants.
Exemple 1 : Optimisation de l'épaisseur de paroi
Ce modèle est un simple caisson surmonté d'un cylindre conique. Les deux parties sont constituées de tôle acier avec différentes épaisseurs de paroi. Celles-ci ont été paramétrées pour réaliser une étude paramétrique. L'objectif est de trouver la configuration dont la contrainte équivalente selon Von Mises ne dépasse pas la limite d'élasticité pour un poids minimal. La charge appliquée est le poids propre ainsi qu'une charge simple composante par composante dans un cas de charge d'exploitation sur une surface rigide au sommet du cône. L'appui est réalisé de manière articulée sur deux lignes à la base du caisson. Le modèle et une vue de l'interface utilisateur avec la paramétrisation peuvent être consultés ci-dessous.
Le script suivant est divisé en sections. Dans la première section, les bibliothèques utilisées sont importées, suivies de la définition des variables. Ici, la limite d'élasticité est indiquée à des fins de visualisation ainsi que la plage de variation des paramètres, qui sont liés à l'épaisseur de tôle. Ensuite, des fonctionnalités sont définies pour la construction de la matrice de variantes et la modification des paramètres globaux dans RFEM 6 dans la suite du processus.
La partie principale du programme commence par la construction de la matrice de variantes des paramètres globaux. Pour cela, toutes les mutations possibles sont générées à partir des variations de paramètres définies au moyen du produit cartésien. Par la suite, celles-ci sont encore limitées selon la condition que l'épaisseur de paroi dans la zone inférieure doit être supérieure ou égale à celle de la partie supérieure.
Ensuite, la connexion avec RFEM 6 est établie via l'API Dlubal. Dans ce cas, vers le modèle actif avec la clé API enregistrée dans le fichier de configuration. Les deux paramètres peuvent également être spécifiés dans l'appel de fonction.
Une boucle est alors exécutée, dans laquelle les paramètres globaux dans RFEM 6 sont ajustés successivement selon la mutation actuelle, le calcul est effectué et les résultats sont lus. Pour forcer une adaptation du modèle et la suppression des résultats, le maillage est supprimé et régénéré avant le calcul. Le maximum de la contrainte équivalente selon Von Mises de toutes les surfaces est déterminé à partir des résultats de la première situation de projet. Le poids total de la structure est déterminé à partir de la force résultante dans la direction z dans le cas de charge de poids propre. À partir de celle-ci, le poids est calculé par division par l'accélération de la pesanteur.
Dans les deux dernières sections, les résultats sont sauvegardés sous forme de tableau Excel et tracés sous forme de diagramme de la contrainte par rapport au poids de la structure.
"""ParaS-VMStress_Thick-min.py
Related RFEM 6 model: "Para_Plate_Joint.rf6"
Model: https://www.dlubal.com/en/downloads-and-information/examples-and-tutorials/models-to-download/006159
Run a parameter study via the Dlubal RFEM API, collect result metrics, export and plot results.
and chart.
Sections:
- Imports
- Configuration
- Functions
- Main execution
- Plotting
"""
# --------------------- Imports ---------------------
import itertools
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from dlubal.api import rfem
# --------------------- Configuration ---------------------
parameter_ranges = {
't_1': {'min': 8, 'max': 12, 'step': 2},
't_2': {'min': 4, 'max': 8, 'step': 2},
} # thicknesses in mm
f_y = 235 # yield strength in MPa
# --------------------- Functions ---------------------
def build_parameter_grid(param_ranges):
"""Build a complete mutation dataset from min/max/step values in mm."""
columns = []
value_lists = []
for name, bounds in param_ranges.items():
min_val = bounds['min']
max_val = bounds['max']
step = bounds['step']
if step <= 0:
raise ValueError(f"Step width for '{name}' must be positive.")
if max_val < min_val:
raise ValueError(f"Max value for '{name}' must be >= min value.")
values = np.arange(min_val, max_val + step * 0.5, step)
columns.append(name)
value_lists.append(values)
combinations = list(itertools.product(*value_lists))
return pd.DataFrame(combinations, columns=columns)
def set_glpa(dl_app, p_name, p_value,):
"""Sets the global parameter p_name to the specified value p_value"""
params = dl_app.get_object_list([rfem.global_parameters.GlobalParameter()])
for p in params:
if p.name == p_name:
p.value = p_value
dl_app.update_object(p)
check=True
return True
if not check:
raise ValueError(f" Error: Parameter '{p_name}' not found!")
# --------------------- Main execution ---------------------
# Build parameter set
para_dfo = build_parameter_grid(parameter_ranges)
para_dfo = para_dfo[para_dfo['t_1'] >= para_dfo['t_2']]
para_df = para_dfo / 1000 # mm -> m
print(f"Parameter set:\n {para_dfo}")
# Connect to RFEM with current model and run parameter study
with rfem.Application() as rf_app:
# Check connection and print model info
app_info = rf_app.get_application_info()
print("Application Info:", app_info)
res_paras=pd.DataFrame(columns=['sigvm','sfz'])
for i in para_df.index:
# set global parameters
print(f"Set mutation {i} with parameters {[x for x in para_df.loc[i]]}")
for p in para_df.columns:
# print(f"Set mutation {i} with parameter {p} to {para_df.loc[i, p]}")
if not set_glpa(rf_app, p, para_df.loc[i, p]):
break
# Remesh
rf_app.delete_mesh()
rf_app.generate_mesh(skip_warnings=True)
# Run calculation
calculation = rf_app.calculate_all(skip_warnings=True)
# Results
if calculation.succeeded:
# Von Mises stresses
res_sigvm = rf_app.get_results(
results_type=rfem.results.STATIC_ANALYSIS_SURFACES_EQUIVALENT_STRESSES_MISES_MESH_NODES,
filters=[
rfem.results.ResultsFilter(column_id='loading', filter_expression='DS1'),
]
).data
sigvm = res_sigvm['sigma_eqv_mises'].max()/10**6 # N/m^2 -> MPa
# self-weight
sfz = rf_app.get_result_table(
table=rfem.results.ResultTable.STATIC_ANALYSIS_SUMMARY_TABLE,
loading=rfem.ObjectId(no=1,object_type=rfem.OBJECT_TYPE_LOAD_CASE)
).data
sfz = float(sfz.loc[6].value) / 10 # N -> kg (g=10 m/s^2 for gravity acceleration)
res_paras.loc[i] = pd.Series({'sigvm':sigvm,'sfz':sfz})
else:
print(f"Calculation failed for mutation {i} with parameters {[x for x in para_df.loc[i]]}")
res_paras.loc[i] = pd.Series({'sigvm':np.nan,'sfz':np.nan})
# Merge parameters and results for export and plotting
para_out = pd.concat([para_dfo,res_paras], axis=1)
# ---------------------- Export Results ---------------------
para_out.to_excel('./ParaS-VMStress_Thick-min_results.xlsx')
# --------------------- Plotting ---------------------
fig, ax = plt.subplots(figsize=(12, 9))
ax.scatter(para_out['sfz'], para_out['sigvm'], color='tab:blue', s=50, label='Mutations')
y_offset = (para_out['sigvm'].max() - para_out['sigvm'].min()) * 0.03
for x, y, mut, t1, t2 in zip(para_out['sfz'], para_out['sigvm'], para_out.index, para_out['t_1'], para_out['t_2']):
ax.text(x , y+y_offset, f"M{int(mut)}, {t1:.0f}, {t2:.0f} mm", fontsize=8, va='top', ha='center')
ax.axhline(y=f_y, color='red', linestyle='--', linewidth=1, label=f'Yield Strength ({f_y} MPa)')
ax.set_title('Stress vs. Mass for Parameter Mutations')
ax.set_xlabel('Mass [kg]')
ax.set_ylabel('Von Mises Stress [MPa]')
ax.legend()
fig.tight_layout()
fig.savefig('./ParaS-VMStress_Thick-min_results.png', dpi=200)
plt.show()Comme le montre cette simple étude paramétrique, par rapport à la configuration initiale avec des épaisseurs de paroi de 12 et 8 mm, des épaisseurs de paroi de 8 mm chacune permettent d'économiser environ 30 % de matériau sans dépasser la limite d'élasticité admissible. Les résultats peuvent être consultés dans l'image suivante.
Exemple 2 : Étude de convergence du maillage
Cet exemple a déjà été utilisé pour l'article technique sur les études de convergence du maillage. Il s'agit ici d'une coque cylindrique pour laquelle une analyse de flambement linéaire est effectuée. Le lien vers l'article technique, qui contient une description plus précise, et le modèle correspondant peuvent être consultés ci-dessous :
L'objectif d'une étude de convergence du maillage est d'adapter la finesse du maillage de manière à ce qu'aucun changement pertinent des résultats ne se produise plus, sans atteindre un nombre d'éléments trop élevé pour permettre un travail économique. Pour cela, un paramètre global (lfe) a été introduit dans le modèle présent. Celui-ci a ensuite été transmis au raffinement du maillage surfacique comme taille de maillage.
Dans le script suivant, l'étude de convergence du maillage est réalisée par modification progressive de la taille de maillage EF souhaitée. Sa structure de base correspond à celle du premier exemple, c'est pourquoi seules les différences sont abordées ici.
Pour analyser l'influence de la finesse du maillage sur le temps de calcul, celui-ci est déterminé à partir du temps de traitement de l'appel de calcul. Les facteurs de charge critiques sont déterminés à partir de l'analyse de stabilité du premier cas de charge. Ensuite, une détermination de la variation du facteur de charge par rapport à l'étape précédente est effectuée afin de pouvoir forcer un arrêt en cas de plage de variation plus grande. Avant la sortie des résultats, une analyse de convergence est effectuée, qui évalue la variation relative des résultats.
"""ParaS-LBA_MeshConvergence-min.py
Related RFEM 6 model: "MeshSensitivity-LBA_AlC.rf6"
Model: https://www.dlubal.com/en/downloads-and-information/examples-and-tutorials/models-to-download/006080
Run a mesh convergence study in RFEM by sweeping the global parameter
`l_fe`, collect critical load factor, compare variants, save and plot results.
Sections:
- Imports
- Configuration
- Main execution
- Plotting
"""
# --------------------- Imports ---------------------
# Standard library
import time
# Third-party
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from dlubal.api import rfem
# ------------------- Configuration -------------------
FE_SIZES_MM = [15, 12, 10, 9, 8, 7, 6, 5, 4, 3] # List of finite element sizes to sweep in millimetres
PARAM_NAME = "l_fe" # Name of the global parameter in the model
CONVERGENCE_THRESHOLD_PCT = 1.0 # Convergence threshold in percentage
FCR_ANA = 1065 # Analytical critical load factor (for reference in the plot)
# ------------------- Functions -------------------
def set_glpa(dl_app, p_name, p_value,):
"""Sets the global parameter p_name to the specified value p_value"""
params = dl_app.get_object_list([rfem.global_parameters.GlobalParameter()])
for p in params:
if p.name == p_name:
p.value = p_value
dl_app.update_object(p)
check=True
return True
if not check:
raise ValueError(f" Error: Parameter '{p_name}' not found!")
# ------------------- Main execution -------------------
with rfem.Application() as rf_app:
# Check connection and print model info
app_info = rf_app.get_application_info()
print("Application Info:", app_info)
results = pd.DataFrame(columns=['l_fe_mm', 'n_elements', 'c_time_s', 'f_cr', 'delta_pct'])
f_prev = None
i = 0
for l_fe in FE_SIZES_MM:
i += 1
l_fe_m = l_fe / 1000.0 # Convert mm to m for use in RFEM
if not set_glpa(dl_app=rf_app, p_name=PARAM_NAME, p_value=l_fe_m):
break
# Remesh
rf_app.delete_mesh()
rf_app.generate_mesh(skip_warnings=True)
n_elem = rf_app.get_mesh_statistics().surface_2D_finite_elements
# Run calculation
t0 = time.perf_counter()
calculation = rf_app.calculate_all(skip_warnings=True)
t1 = time.perf_counter()
c_time = t1 - t0
# Results
if calculation.succeeded:
f_cr = rf_app.get_results(
results_type=rfem.results.STABILITY_ANALYSIS_CRITICAL_LOAD_FACTORS,
filters=[
rfem.results.ResultsFilter(column_id='loading', filter_expression='LC1'),
]
).data.loc[0].f
if f_prev is not None and f_prev != 0 and f_cr is not None:
delta_pct = abs(f_cr - f_prev) / abs(f_prev) * 100.0
delta_str = f"{delta_pct:>4.2f}"
else:
delta_pct = None
delta_str = f"{'---':>4}"
if f_cr is not None:
print(f"Mutation {i}: l_fe = {l_fe:.1f} mm | n_elements = {n_elem} | f_cr = {f_cr:.2f} | delta = {delta_str} %")
results.loc[i] = pd.Series({
'l_fe_mm': l_fe, 'n_elements': n_elem, 'c_time_s': c_time,
'f_cr': f_cr, 'delta_pct': delta_pct
})
f_prev = f_cr
else:
print(f"Calculation failed for mutation {i} with mesh size {l_fe:.1f} mm")
results.loc[i] = pd.Series({
'l_fe_mm': l_fe, 'n_elements': n_elem, 'c_time_s': np.nan,
'f_cr': np.nan, 'delta_pct': np.nan
})
# Calculate 1 / n_elements for the second plot
results['inv_n_elements'] = 1.0 / results['n_elements']
# Calculate relative deviation from minimum critical load factor
results['rel_dev_f_cr_min'] = (results['f_cr'] - results['f_cr'].min()) / results['f_cr'].min() * 100
# Convergence analysis
print(f"\n Convergence analysis (threshold: < {CONVERGENCE_THRESHOLD_PCT} % change):")
convergence_point = None
for (l_mm, n, a, d) in results.loc[:, ['l_fe_mm', 'n_elements', 'f_cr', 'delta_pct']].itertuples(index=False):
if d is not None and d < CONVERGENCE_THRESHOLD_PCT:
convergence_point = (l_mm, n, a)
break
if convergence_point:
print(f" Convergence reached at l_fe = {convergence_point[0]:.1f} mm")
print(f" f_cr = {convergence_point[2]:.2f} | n_elements = {convergence_point[1]}")
else:
print(" Convergence not reached yet – finer mesh required.")
# ---------------------- Export Results ---------------------
results.to_excel('./ParaS-LBA_MeshConvergence-min_results.xlsx')
# ------------------- Plotting -------------------
fig, axes = plt.subplots(2, 2, figsize=(14, 10))
# Plot 1: f_cr vs. n_elements
ax1 = axes[0, 0]
ax1.plot(results['n_elements'], results['f_cr'], 'o-', color='tab:blue', markersize=6)
ax1.axhline(y=FCR_ANA, color='red', linestyle='--', linewidth=1, label=f'$F_{{cr,ana}}$ = {FCR_ANA} kN')
ax1.set_xlabel('Number of Elements')
ax1.set_ylabel('Critical Load Factor [kN]')
ax1.set_title('Critical Load Factor vs. Number of Elements')
ax1.grid(True, alpha=0.3)
ax1.legend()
# Plot 2: f_cr vs. 1/n_elements
ax2 = axes[0, 1]
ax2.plot(results['inv_n_elements'], results['f_cr'], 'o-', color='tab:blue', markersize=6)
ax2.axhline(y=FCR_ANA, color='red', linestyle='--', linewidth=1, label=f'$F_{{cr,ana}}$ = {FCR_ANA} kN')
ax2.set_xlabel('1 / Number of Elements')
ax2.set_ylabel('Critical Load Factor [kN]')
ax2.set_title('Critical Load Factor vs. 1 / Number of Elements')
ax2.grid(True, alpha=0.3)
ax2.legend()
# Plot 3: Relative Deviation vs. Calculation Time
ax3 = axes[1, 0]
ax3.plot(results['c_time_s'], results['rel_dev_f_cr_min'], 'o-', color='tab:blue', markersize=6)
ax3.set_xlabel('Calculation Time [s]')
ax3.set_ylabel('Relative Deviation [%]')
ax3.set_title('Relative Deviation vs. Calculation Time')
ax3.grid(True, alpha=0.3)
# Plot 4: Relative Deviation vs. l_fe_mm
ax4 = axes[1, 1]
ax4.plot(results['l_fe_mm'], results['rel_dev_f_cr_min'], 'o-', color='tab:blue', markersize=6)
ax4.set_xlabel('FE Mesh Size [mm]')
ax4.set_ylabel('Relative Deviation [%]')
ax4.set_title('Relative Deviation vs. FE Mesh Size')
ax4.grid(True, alpha=0.3)
fig.suptitle('Mesh Convergence Analysis', fontsize=14, fontweight='bold')
fig.tight_layout()
fig.savefig('./ParaS-LBA_MeshConvergence-min_results.png', dpi=200)
plt.show()L'image suivante montre les diagrammes générés par le script du programme, tels qu'ils sont également présentés et évalués plus en détail dans l'article technique mentionné au début.
