.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "_gallery/AxialStretchingCase/run_axial_stretching.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__gallery_AxialStretchingCase_run_axial_stretching.py: Axial Stretching ================ This case tests the axial stretching of a rod. The expected behavior is supposed to be like a spring-gravity motion, but with a rod. A rod is fixed at one end and a force is applied at the other end. The rod stretches and the displacement of the tip is compared with the analytical solution. .. GENERATED FROM PYTHON SOURCE LINES 11-19 .. code-block:: Python # isort:skip_file import numpy as np from matplotlib import pyplot as plt import elastica as ea .. GENERATED FROM PYTHON SOURCE LINES 20-24 Simulation Setup ---------------- We define a simulator class that inherits from the necessary mixins. This makes constraints, forces, and damping evailable to the system. .. GENERATED FROM PYTHON SOURCE LINES 24-35 .. code-block:: Python class StretchingBeamSimulator( ea.BaseSystemCollection, ea.Constraints, ea.Forcing, ea.Damping, ea.CallBacks ): pass stretch_sim = StretchingBeamSimulator() final_time = 200.0 .. GENERATED FROM PYTHON SOURCE LINES 36-43 Rod Setup --------- Next, we set up the test parameters for the simulating rods. This includes the number of elements, the start position, direction, normal, length, radius, density, and Young's modulus of the rod. For this case, we have fixed boundary condition at one end, and we apply external force at the other end. .. GENERATED FROM PYTHON SOURCE LINES 43-82 .. code-block:: Python # setting up test params n_elem = 19 start = np.zeros((3,)) direction = np.array([1.0, 0.0, 0.0]) normal = np.array([0.0, 1.0, 0.0]) base_length = 1.0 base_radius = 0.025 base_area = np.pi * base_radius**2 density = 1000 youngs_modulus = 1e4 # For shear modulus of 1e4, nu is 99! poisson_ratio = 0.5 shear_modulus = youngs_modulus / (poisson_ratio + 1.0) stretchable_rod = ea.CosseratRod.straight_rod( n_elem, start, direction, normal, base_length, base_radius, density, youngs_modulus=youngs_modulus, shear_modulus=shear_modulus, ) stretch_sim.append(stretchable_rod) stretch_sim.constrain(stretchable_rod).using( ea.OneEndFixedBC, constrained_position_idx=(0,), constrained_director_idx=(0,) ) end_force_x = 1.0 end_force = np.array([end_force_x, 0.0, 0.0]) stretch_sim.add_forcing_to(stretchable_rod).using( ea.EndpointForces, 0.0 * end_force, end_force, ramp_up_time=1e-2 ) .. GENERATED FROM PYTHON SOURCE LINES 83-85 Damping is added to the system to help it reach a steady state. We use an `AnalyticalLinearDamper` to add damping to the rod. .. GENERATED FROM PYTHON SOURCE LINES 85-97 .. code-block:: Python # add damping dl = base_length / n_elem dt = 0.1 * dl damping_constant = 0.1 stretch_sim.dampen(stretchable_rod).using( ea.AnalyticalLinearDamper, damping_constant=damping_constant, time_step=dt, ) .. GENERATED FROM PYTHON SOURCE LINES 98-102 Callbacks --------- A callback object is passed to the simulator to record states of the rod during the simulation. This is useful for post-processing the results. .. GENERATED FROM PYTHON SOURCE LINES 102-137 .. code-block:: Python # Add call backs class AxialStretchingCallBack(ea.CallBackBaseClass): """ Tracks the velocity norms of the rod """ def __init__(self, step_skip: int, callback_params: dict) -> None: super().__init__() self.every = step_skip self.callback_params = callback_params def make_callback( self, system: ea.typing.RodType, time: np.float64, current_step: int ) -> None: if current_step % self.every == 0: self.callback_params["time"].append(time) # Collect only x self.callback_params["position"].append( system.position_collection[0, -1].copy() ) self.callback_params["velocity_norms"].append( np.linalg.norm(system.velocity_collection.copy()) ) return recorded_history: dict[str, list] = ea.defaultdict(list) stretch_sim.collect_diagnostics(stretchable_rod).using( AxialStretchingCallBack, step_skip=200, callback_params=recorded_history ) .. GENERATED FROM PYTHON SOURCE LINES 138-142 Finalize and Run ---------------- We finalize the simulator and create the time-stepper. The `PositionVerlet` time-stepper is used to integrate the system. .. GENERATED FROM PYTHON SOURCE LINES 142-154 .. code-block:: Python stretch_sim.finalize() timestepper: ea.typing.StepperProtocol = ea.PositionVerlet() # timestepper = PEFRL() total_steps = int(final_time / dt) print("Total steps", total_steps) dt = final_time / total_steps time = 0.0 for i in range(total_steps): time = timestepper.step(stretch_sim, time, dt) .. rst-class:: sphx-glr-script-out .. code-block:: none Total steps 38000 .. GENERATED FROM PYTHON SOURCE LINES 155-160 Post-Processing --------------- Finally, we plot the results and compare them with the analytical solution. The analytical solution is calculated using the first-order theory with both the base length and the modified length. .. GENERATED FROM PYTHON SOURCE LINES 160-176 .. code-block:: Python # First-order theory with base-length expected_tip_disp = end_force_x * base_length / base_area / youngs_modulus # First-order theory with modified-length, gives better estimates expected_tip_disp_improved = ( end_force_x * base_length / (base_area * youngs_modulus - end_force_x) ) fig = plt.figure(figsize=(10, 8), frameon=True, dpi=150) ax = fig.add_subplot(111) ax.plot(recorded_history["time"], recorded_history["position"], lw=2.0) ax.hlines(base_length + expected_tip_disp, 0.0, final_time, "k", "dashdot", lw=1.0) ax.hlines( base_length + expected_tip_disp_improved, 0.0, final_time, "k", "dashed", lw=2.0 ) plt.show() .. image-sg:: /_gallery/AxialStretchingCase/images/sphx_glr_run_axial_stretching_001.png :alt: run axial stretching :srcset: /_gallery/AxialStretchingCase/images/sphx_glr_run_axial_stretching_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 2.276 seconds) .. _sphx_glr_download__gallery_AxialStretchingCase_run_axial_stretching.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: run_axial_stretching.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: run_axial_stretching.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: run_axial_stretching.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_