.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "advanced_examples/304L_viscoplastic_calibration/plot_304L_a_data_analysis.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_advanced_examples_304L_viscoplastic_calibration_plot_304L_a_data_analysis.py: 304L bar data analysis ---------------------- In this example, we analyze the tension data for 304L at several strain rates from :cite:p:`laser_weld_paper`. Since 304L material is known to exhibit rate dependence that is measurable when increasing strain rates an order of magnitude or more, this data set is a good data set to analyze for rate dependence. After analyzing the data, we will then use MatCal tools to help chose a model for for the rate dependence. .. note:: Useful Documentation links: #. :class:`~matcal.core.data_importer.BatchDataImporter` #. :meth:`~matcal.core.data.DataCollection.plot` #. :func:`~matcal.core.data_analysis.determine_pt2_offset_yield` First we import the python libraries we will use. For this example, we will use MatCal, NumPy and Matplotlib. .. GENERATED FROM PYTHON SOURCE LINES 21-30 .. code-block:: Python from matcal import * import matplotlib.pyplot as plt import numpy as np plt.rc('text', usetex=True) plt.rc('font', family='serif') plt.rc('font', size=12) figsize = (4,3) .. GENERATED FROM PYTHON SOURCE LINES 31-35 Next, we import the data using the :class:`~matcal.core.data_importer.BatchDataImporter`. These files have been preprocessed to include state information in each file, so that they can be imported with `engineering_strain_rate` and `temperature` state variables pre-assigned to each experiment. .. GENERATED FROM PYTHON SOURCE LINES 35-37 .. code-block:: Python data = BatchDataImporter("ductile_failure_small_tension_data/*.csv").batch .. GENERATED FROM PYTHON SOURCE LINES 38-49 Now that the data is imported, we create data structures that will aid in plotting the yield stresses and stresses at 5 percent strain as a function of rate. We look at these stresses to determine if both yield and hardening are rate dependent. So for each experiment in the data set, we save these values in our data structures. MatCal's :func:`~matcal.core.data_analysis.determine_pt2_offset_yield` function is useful here for extracting the 0.2% offset yield stress from engineering stress strain curves. .. GENERATED FROM PYTHON SOURCE LINES 49-69 .. code-block:: Python yield_stresses = {0.0001:[], 0.01:[], 500.0:[], 1800.0:[], 3600.0:[]} five_percent_strain_stresses = {0.0001:[], 0.01:[], 500.0:[], 1800.0:[], 3600.0:[]} for state, data_sets in data.items(): for data in data_sets: rate = state["engineering_strain_rate"] yield_pt = determine_pt2_offset_yield(data, 29e3) yield_stresses[rate].append(yield_pt[1]) five_strain = np.interp(0.05, data["engineering_strain"], data["engineering_stress"]) five_percent_strain_stresses[rate].append(five_strain) .. GENERATED FROM PYTHON SOURCE LINES 70-74 With the data organized as needed, we can create box blots of these values at each rate. This will allow us to see how these values change for each of the measured rates. .. GENERATED FROM PYTHON SOURCE LINES 74-84 .. code-block:: Python plt.figure(figsize=figsize, constrained_layout=True) bp_yield=plt.boxplot(yield_stresses.values(), labels=yield_stresses.keys()) plt.setp(bp_yield['boxes'], color='tab:blue') plt.xlabel("engineering strain rate (1/s)") plt.ylabel("stress (ksi)") bp_5 =plt.boxplot(five_percent_strain_stresses.values(), labels=five_percent_strain_stresses.keys()) plt.setp(bp_5['boxes'], color='tab:red') plt.legend([bp_yield["boxes"][0], bp_5["boxes"][0]], ['yield stress', 'stress at 5\% strain']) .. image-sg:: /advanced_examples/304L_viscoplastic_calibration/images/sphx_glr_plot_304L_a_data_analysis_001.png :alt: plot 304L a data analysis :srcset: /advanced_examples/304L_viscoplastic_calibration/images/sphx_glr_plot_304L_a_data_analysis_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none /gpfs/knkarls/projects/matcal-stable/external_matcal/documentation/advanced_examples/304L_viscoplastic_calibration/plot_304L_a_data_analysis.py:83: SyntaxWarning: invalid escape sequence '\%' 'stress at 5\% strain']) .. GENERATED FROM PYTHON SOURCE LINES 85-124 From these plots, we can see that the material does exhibit rate dependence when the engineering strain rate changes several orders of magnitude. Rate dependence in the material yield is clear. The stresses at five percent strain show that the material hardening is likely rate independent. This is apparent because the stress increase at the different rates does not increase at the stress at 5% strain as would be expected if the material hardening was also rate dependent. Instead it decreases. This decrease is likely due to heating due to plastic work in the material at high rate. We now plot just the yield data on a ``semilogx`` plot to visually assess the relationship between yield stress and strain rate. Two commonly used strain rate dependent yield models for metals include the Johnson-Cook model (JC) and the Power-law Breakdown model (PLB). The functional form for the JC rate dependence model is: .. math:: Y\left(\dot{\epsilon}^p\right) = Y_0\left[1+C\ln\left(\frac{\dot{\epsilon}^p} {\dot{\epsilon}_0}\right)\right] where :math:`Y_0` is the rate independent yield stress, :math:`C` is a calibration constant, :math:`\dot{\epsilon}^p` is the material plastic strain rate, and :math:`\dot{\epsilon}_0` is a reference strain rate under which the material is rate independent. The functional form for the PLB rate dependence model is: .. math:: Y\left(\dot{\epsilon}^p\right) = Y_0\left[1+\text{sinh}^{-1}\left(\left( \frac{\dot{\epsilon}^p}{g}\right)^{(1/m)}\right)\right] where :math:`Y_0` is the rate independent yield stress, and :math:`g` and :math:`b` are a calibration constants. As a function of strain rate, the JC model is linear in a ``semilogx`` plot while the PLB can exhibit curvature. To see how these data looks and if one of these model are clearly more appropriate, we plot the data on a ``semilogx`` First we put the data into a MatCal :class:`matcal.core.data.Data` class. .. GENERATED FROM PYTHON SOURCE LINES 124-131 .. code-block:: Python yield_dc = DataCollection("yeild vs rate") for rate in yield_stresses: rate_state = State(f"rate_{rate}", rate=rate) for yield_stress in yield_stresses[rate]: data = convert_dictionary_to_data({"yield":[yield_stress]}) data.set_state(rate_state) yield_dc.add(data) .. GENERATED FROM PYTHON SOURCE LINES 132-133 Next, we plot the data. .. GENERATED FROM PYTHON SOURCE LINES 133-154 .. code-block:: Python plt.figure(figsize=figsize, constrained_layout=True) def plot_dc_by_state(data_collection, label=None, color=None, best_index=None): for state in data_collection: if best_index is None: for idx, data in enumerate(data_collection[state]): plt.semilogx(state["rate"], data["yield"][0], marker='o', label=label, color=color) if color is not None and label is not None: label = "_"+label else: data = data_collection[state][best_index] plt.semilogx(state["rate"], data["yield"][0], marker='o', label=label, color=color) if color is not None and label is not None: label = "_"+label plt.xlabel("engineering strain rate (1/s)") plt.ylabel("yield stress (ksi)") plot_dc_by_state(yield_dc) plt.show() .. image-sg:: /advanced_examples/304L_viscoplastic_calibration/images/sphx_glr_plot_304L_a_data_analysis_002.png :alt: plot 304L a data analysis :srcset: /advanced_examples/304L_viscoplastic_calibration/images/sphx_glr_plot_304L_a_data_analysis_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 155-170 Upon inspection, it is not immediately clear which model will fit the data better. This is likely due to the fact that there is significant scatter at the different strain rates and there is no data in the intermediate strain rates. As a result, we will use MatCal tools to help decide which model we should choose based on these data. To begin, we calibrate each python model to these data. We already have our data, so we need to create models that can predict the trend in the data. We use MatCal's :class:`~matcal.core.models.PythonModel` to implement the models using python functions. These two models are defined below. First, we define the JC model python function. .. GENERATED FROM PYTHON SOURCE LINES 170-175 .. code-block:: Python def jc_rate_dependence_model(Y_0, C, ref_strain_rate, rate): yield_stresses = np.atleast_1d(Y_0*(1+C*np.log(rate/ref_strain_rate))) yield_stresses[np.atleast_1d(rate) < ref_strain_rate] = Y_0 return {"yield":yield_stresses} .. GENERATED FROM PYTHON SOURCE LINES 176-179 Then we create the python model, name it and add the reference strain rate in as a state parameter as it will be uncalibrated. .. GENERATED FROM PYTHON SOURCE LINES 179-183 .. code-block:: Python jc_rate_model = PythonModel(jc_rate_dependence_model) jc_rate_model.set_name("python_jc_rate_model") jc_rate_model.add_constants(ref_strain_rate=1e-5) .. GENERATED FROM PYTHON SOURCE LINES 184-185 Next we define the PLB model. .. GENERATED FROM PYTHON SOURCE LINES 185-189 .. code-block:: Python def plb_rate_dependence_model(Y_0, g_star, m, rate): yield_stress = Y_0*(1+np.arcsinh((rate/10**(g_star))**(1/m))) return {"yield":np.atleast_1d(yield_stress)} .. GENERATED FROM PYTHON SOURCE LINES 190-196 Note that we are calibrating the constant :math:`g` on a log scale. We create a parameter :math:`g^*`, such that :math:`g=10^{g*}`. This is needed because the model would otherwise appear insensitive to :math:`g` in MatCal studies. With the function created, we now make the model and name it. .. GENERATED FROM PYTHON SOURCE LINES 196-199 .. code-block:: Python plb_rate_model = PythonModel(plb_rate_dependence_model) plb_rate_model.set_name("python_plb_rate_model") .. GENERATED FROM PYTHON SOURCE LINES 200-204 Next, we now define the parameters that will be calibrated for these model. Both models will need to calibrate the rate independent yield stress, :math:`Y_0`. The bounds for this parameter are set based on what we observe in the low strain rate data. .. GENERATED FROM PYTHON SOURCE LINES 204-206 .. code-block:: Python Y_0 = Parameter("Y_0", 20, 60) .. GENERATED FROM PYTHON SOURCE LINES 207-216 The PLB model requires the two parameters :math:`g_star` and :math:`m`. The :math:`g` parameter is meant to represent a change in behavior for the rate dependence. It is based on experimental observations that show materials sensitivity to rate increase at high strain rates. In this model, it can be considered a reference rate above which the material becomes more rate sensitive. Since this is generally at higher rates for metals, we restrict this reference rate to be between 100 and 10000 per second. Any lower or higher, and the parameter will be used in an unintended fashion. The bounds for :math:`m` are set based on our previous experience with the model for austenitic stainless steels. .. GENERATED FROM PYTHON SOURCE LINES 216-219 .. code-block:: Python g_star = Parameter("g_star", 2, 4) m = Parameter("m", 2, 15) .. GENERATED FROM PYTHON SOURCE LINES 220-223 The only unique parameter for the JC model is the calibration parameter :math:`C` which we also set bounds for using our previous experience with the model for metals. .. GENERATED FROM PYTHON SOURCE LINES 223-225 .. code-block:: Python C = Parameter("C", 0.001, 0.1) .. GENERATED FROM PYTHON SOURCE LINES 226-229 The last component needed for our calibrations is the objective. We use a :class:`~matcal.core.objective.Objective` to fit our python models to the data. .. GENERATED FROM PYTHON SOURCE LINES 229-231 .. code-block:: Python obj = Objective("yield") .. GENERATED FROM PYTHON SOURCE LINES 232-233 Now we can setup our calibrations and save the results. .. GENERATED FROM PYTHON SOURCE LINES 233-243 .. code-block:: Python jc_cal = GradientCalibrationStudy(Y_0, C) jc_cal.add_evaluation_set(jc_rate_model, obj, yield_dc) jc_cal.set_working_directory("jc", remove_existing=True) jc_cal_results = jc_cal.launch() plb_cal = GradientCalibrationStudy(Y_0, g_star, m) plb_cal.add_evaluation_set(plb_rate_model, obj, yield_dc) plb_cal.set_working_directory("plb", remove_existing=True) plb_cal_results = plb_cal.launch() .. GENERATED FROM PYTHON SOURCE LINES 244-246 With the calibrations complete, we plot the fits against the data and print the best fit information. .. GENERATED FROM PYTHON SOURCE LINES 246-265 .. code-block:: Python jc_best_idx = jc_cal_results.best_evaluation_index jc_best_sim = jc_cal_results.simulation_history[jc_rate_model.name] plb_best_idx = plb_cal_results.best_evaluation_index plb_best_sim = plb_cal_results.simulation_history[plb_rate_model.name] plt.figure(figsize=figsize, constrained_layout=True) plot_dc_by_state(yield_dc, "experiments" , 'k') rates=np.logspace(-4,4,100) plt.semilogx(rates, jc_rate_dependence_model(**jc_cal_results.best.to_dict(), ref_strain_rate=1e-5, rate=rates)['yield'], label="JC model") plot_dc_by_state(jc_best_sim, None, 'tab:blue', best_index=jc_best_idx) plt.semilogx(rates, plb_rate_dependence_model(**plb_cal_results.best.to_dict(), rate=rates)['yield'], label="PLB model") plot_dc_by_state(plb_best_sim, None, 'tab:orange', best_index=plb_best_idx) plt.legend() plt.show() .. image-sg:: /advanced_examples/304L_viscoplastic_calibration/images/sphx_glr_plot_304L_a_data_analysis_003.png :alt: plot 304L a data analysis :srcset: /advanced_examples/304L_viscoplastic_calibration/images/sphx_glr_plot_304L_a_data_analysis_003.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 266-270 From this plot, both models appear to represent the data well. Looking at the final total objective and parameters will reveal which model fits the data better and if there were any issues in the fitting process. .. GENERATED FROM PYTHON SOURCE LINES 270-280 .. code-block:: Python jc_eval_set_name = f"{jc_rate_model.name}:{obj.name}" jc_best_obj = jc_cal_results.best_total_objective print("JC model best fit objective:", jc_best_obj) print(jc_cal_results.best.to_dict(),"\n") plb_eval_set_name = f"{plb_rate_model.name}:{obj.name}" plb_best_obj = plb_cal_results.best_total_objective print("PLB model best fit objective:", plb_best_obj) print(plb_cal_results.best.to_dict(),"\n") .. rst-class:: sphx-glr-script-out .. code-block:: none JC model best fit objective: 0.004652787936813707 OrderedDict({'Y_0': 32.444999516, 'C': 0.039894200734}) PLB model best fit objective: 0.004014596531156934 OrderedDict({'Y_0': 33.037043219, 'g_star': 4.0, 'm': 8.3580048779}) .. GENERATED FROM PYTHON SOURCE LINES 281-295 The objectives show that the PLB model provides a better fit to these data. However, we can see in the fit that its parameter for :math:`g_star` is hitting its upper bound. This is showing that the calibration is adjusting that parameter outside is intended use case for this model and is the first indication that we should use JC over PLB. To look into this further, we will perform a sensitivity study on the objective with respect to the parameters in both models. If the objective is not very sensitive to the parameters, it can be an indication that the parameter is not very well defined by the given objective. We start with the JC model. We need to redefine our parameters and add distributions to them to support the sensitivity study. We will assign them a ``uniform_uncertain`` distribution. .. GENERATED FROM PYTHON SOURCE LINES 295-298 .. code-block:: Python Y_0 = Parameter("Y_0", 20, 60, distribution="uniform_uncertain") C = Parameter("C", 0.001, 0.1, distribution="uniform_uncertain") .. GENERATED FROM PYTHON SOURCE LINES 299-309 Now we can create our :class:`~matcal.dakota.sensitivity_studies.LhsSensitivityStudy` and add our evaluation set. .. note:: We update the metric function of the objective to be the :class:`~matcal.core.objective.L2NormMetricFunction`. This can provide a more interpretable result for sensitivity analyses than the :class:`~matcal.core.objective.SumSquaresMetricFunction`. .. GENERATED FROM PYTHON SOURCE LINES 309-314 .. code-block:: Python sens = LhsSensitivityStudy(Y_0, C) obj.set_metric_function(L2NormMetricFunction()) sens.add_evaluation_set(jc_rate_model, obj, yield_dc) sens.set_working_directory("jc_sens", remove_existing=True) .. GENERATED FROM PYTHON SOURCE LINES 315-319 We want to perform the study on the overall objective and to perform a Sobol index study. The Sobol index study provides a global sensitivity of the models to the input parameters and are valid for nonlinear responses. .. GENERATED FROM PYTHON SOURCE LINES 319-322 .. code-block:: Python sens.use_overall_objective() sens.make_sobol_index_study() .. GENERATED FROM PYTHON SOURCE LINES 323-327 Since the python model is inexpensive, we take 2500 samples. For a problem with a more computationally expensive model, you should run many studies with increasing samples until the Sobol indices converge. .. GENERATED FROM PYTHON SOURCE LINES 327-330 .. code-block:: Python num_samples = 2500 sens.set_number_of_samples(num_samples) .. GENERATED FROM PYTHON SOURCE LINES 331-335 For python models, performance gains can be achieved by running the evaluations in serial. Parallel evaluations require additional overhead that can decrease performance when using inexpensive models such as a :class:`~matcal.core.models.PythonModel`. .. GENERATED FROM PYTHON SOURCE LINES 335-337 .. code-block:: Python sens.run_in_serial() .. GENERATED FROM PYTHON SOURCE LINES 338-339 We now launch the study and save the results. .. GENERATED FROM PYTHON SOURCE LINES 339-341 .. code-block:: Python jc_sens_results = sens.launch() .. GENERATED FROM PYTHON SOURCE LINES 342-343 The sensitivity study for the PLB model is setup the same way. .. GENERATED FROM PYTHON SOURCE LINES 343-354 .. code-block:: Python g_star = Parameter("g_star", 2, 4,distribution="uniform_uncertain") m = Parameter("m", 2, 15, distribution="uniform_uncertain") sens = LhsSensitivityStudy(Y_0, g_star, m) sens.add_evaluation_set(plb_rate_model, obj, yield_dc) sens.set_number_of_samples(num_samples) sens.set_working_directory("plb_sens", remove_existing=True) sens.use_overall_objective() sens.make_sobol_index_study() sens.run_in_serial() .. GENERATED FROM PYTHON SOURCE LINES 355-356 We once again launch the study and save the results. .. GENERATED FROM PYTHON SOURCE LINES 356-358 .. code-block:: Python plb_sens_results = sens.launch() .. GENERATED FROM PYTHON SOURCE LINES 359-360 With both studies complete, we can print the results for analysis. .. GENERATED FROM PYTHON SOURCE LINES 360-363 .. code-block:: Python print("JC sensitivity results:", jc_sens_results.sobol) print("PLB sensitivity results:", plb_sens_results.sobol) .. rst-class:: sphx-glr-script-out .. code-block:: none JC sensitivity results: Y_0: [[0.48219582 0.7043585 ]] C: [[0.25280175 0.49818737]] PLB sensitivity results: Y_0: [[0.85005052 0.91571222]] g_star: [[0.02345957 0.10134016]] m: [[0.00910786 0.04080117]] .. GENERATED FROM PYTHON SOURCE LINES 364-380 The results are printed for each parameter and show the main effects and total effects, respectively. These responses show that the JC model should be used since both its parameters have significant main and total effects. With main effects of ~> 0.3 and total effects ~> 0.6, Both :math:`Y_0` and :math:`C` have a measurable influence on the objective and are significantly coupled. While for the PLB model, :math:`m` and :math:`g_star` only have a slight influence on the objective with main effects < 0.1. Even the total effect of :math:`m` is < 0.1. This indicates they cannot be well calibrated over their expected ranges with the available data and this model should not be used. To finish this study, we save the JC calibrated parameters and the yield vs rate data to a file for use in the full finite element calibration. .. GENERATED FROM PYTHON SOURCE LINES 380-381 .. code-block:: Python matcal_save("JC_parameters.serialized", jc_cal_results.best.to_dict()) matcal_save("rate_data.joblib", yield_dc) .. rst-class:: sphx-glr-timing **Total running time of the script:** (35 minutes 22.710 seconds) .. _sphx_glr_download_advanced_examples_304L_viscoplastic_calibration_plot_304L_a_data_analysis.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_304L_a_data_analysis.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_304L_a_data_analysis.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_304L_a_data_analysis.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_