Search code examples
pythonabaqus

Part (or instance) rotation in abaqus with Euler angle

I tried to model composite with many inclusions (cylinder). I need to rotate the cylinders to the proper position with given Euler angles. But I found from abaqus documentation that instance in abaqus can only be rotated by certain axis.

So I was wondering is there a way in abaqus (or abaqus-python) to rotate instance by Euler angles?

Solution

  • Answered for Abaqus v6.14

    You can do this using Abaqus Python commands in CAE.

    Recall the transformation matrices for some rotations of an angle θ of objects about axes x, y, and z:

    enter image description here

    A rotation about some arbitrary axis u looks like

    enter image description here

    Euler angle transformations look something like the product

    enter image description here

    or similar, depending on whether you follow Proper Euler angle conventions or Tait-Bryan conventions. The above is a Tait-Bryan sequence which rotates an object about the x axis by angle γ, then y axis about angle β, then z axis about angle α.

    Abaqus rotations are defined in CAE via the command ra.rotate(instanceList, axisPoint, axisDirection, angle), where

    • ra is an assembly object (mdb.models[...].rootAssembly)
    • instanceList is an iterable* (tuple, list, or similar) of strings, specifying the instances to be rotated,
    • axisPoint is an iterable of floats specifying the starting point of the axis to be rotated around (x,y,z),
    • axisDirection is an iterable of floats specifying the direction of the axis,
    • angle is the (positive anticlockwise) angle you wish to rotate the object by, in degrees.

    For example, if I wanted to rotate an instance named** with the string 'Part-1-MESH-1-1' under some assembly ra about an axis with origin (a, b, c) in the direction (u, v, w) about an anticlockwise angle t degrees, I would run the command

    ra.rotate(instanceList=('PART-1-MESH-1-1', ), axisPoint=(a, b, c), axisDirection=(u, v, w), angle=t)

    This means that in CAE, the Euler angle transformations can be done sequentially in 3 rotate commands.

    Your script (incl. helper functions) should look something like this:

    """
Assumes a Right-handed coordinate system.
"""

import numpy as np


def _Ru(t, u):
    return np.array([[np.cos(t) + (u[0]**2)*(1 - np.cos(t)), u[0]*u[1]*(1 - np.cos(t)) - u[2]*np.sin(t), u[0]*u[2]*(1 - np.cos(t)) + u[1]*np.sin(t)],
                     [u[1]*u[0]*(1 - np.cos(t)) + u[2]*np.sin(t), np.cos(t) + (u[1]**2)*(1 - np.cos(t)), u[1]*u[2]*(1 - np.cos(t)) - u[0]*np.sin(t)],
                     [u[2]*u[0]*(1 - np.cos(t)) - u[1]*np.sin(t), u[2]*u[1]*(1 - np.cos(t)) + u[0]*np.sin(t), np.cos(t) + (u[2]**2)*(1 - np.cos(t))]])

def general_rotation(axis_order, angles, axes_dir=((1.0,0.0,0.0),(0.0,1.0,0.0))):
    """
    Inputs:
        axis_order - Iterable of ints specifying the axis order to apply the transformation.
                     For example, (0,1,2) means to apply the rotations around the (x,y,z) axes, in that order.
        angles     - Iterable of angles (in degrees) specifying the rotations. This should be the same length as axis_order.
        axes_dir   - Iterable of iterable of floats, specifying two orthogonal directions forming the local x and y axes, respectively. Defaults to global x and y axes.
    """

    # Convert all angles to radians
    angles = np.deg2rad(angles)

    # Calculate the third (z) axis and normalise all axes
    ax0 = np.array(axes_dir[0])
    ax1 = np.array(axes_dir[1])
    ax2 = np.cross(ax0, ax1)
    ax0 = ax0/np.linalg.norm(ax0)
    ax1 = ax1/np.linalg.norm(ax1)
    ax2 = ax2/np.linalg.norm(ax2)
    ax = [ax0, ax1, ax2] # Or similar iterable but must be mutable

    intermediate_axes = [None]*len(axis_order)

    # Calculate total transformation
    for i, a in enumerate(axis_order):

        # Store intermediate axes
        intermediate_axes[i] = list(ax)

        R = _Ru(angles[i], ax[a])

        # Don't bother transforming current axis that's being rotated about
        zot = [0,1,2]
        zot.remove(a)
        for aj in zot:
            ax[aj] = np.dot(R, ax[aj])

    return intermediate_axes

    You can get the output and rotate the instance via

    rotate_axes = general_rotation(axis_order, angles, axes_dir)
for i in range(len(rotate_axes)):
    ra.rotate(instanceList=(instanceName, ), axisPoint=origin, axisDirection=rotate_axes[i][axis_order[i]], angle=angles[i])

    With ra as the assembly object, instanceName as a string representing the name of the instance to be rotated, and origin as the origin of the axes.

    *In Python, when specifying tuple objects with only one member, a comma must be put after the only member. E.g. a tuple containing the number 1 is constructed by (1,), not (1). You can also use lists instead.

    **The exact string can be found in the model tree:

    enter image description here