connectorOutputVariablesTest.py

You can view and download this file on Github: connectorOutputVariablesTest.py

  1#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  2# This is an EXUDYN example
  3#
  4# Details:  Reading the output variables of connectors and joints (#2735), and what happens with one an
  5#           item does not declare:
  6#           (1) ObjectContactCoordinate: a 1D mass resting on a stop; Distance is the gap, negative in
  7#               contact, and the contact force holds the mass; with activeConnector=False the mass falls
  8#               through the stop;
  9#           (2) ObjectJointRevolute2D: the Displacement between the two joint points stays zero;
 10#           (3) ObjectJointPrismatic2D with the rotation constrained: Distance along the axis (the axis
 11#               need not be normalized) and Rotation, which stays zero;
 12#           (4) ObjectJointPrismatic2D with the rotation free: Rotation is the rotation of the body
 13#               relative to the ground;
 14#           (5) an output variable that an item does not declare raises an exception - Rotation of
 15#               ObjectJointRevolute2D, and ObjectContactCircleCable2D, which declares none.
 16#
 17# Author:   Johannes Gerstmayr
 18# Date:     2026-10-01
 19#
 20# Copyright:This file is part of Exudyn. Exudyn is free software. You can redistribute it and/or modify it under the terms of the Exudyn license. See 'LICENSE.txt' for more details.
 21#
 22#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 23
 24import exudyn as exu
 25from exudyn.utilities import * #includes itemInterface and rigidBodyUtilities
 26from exudyn.beams import GenerateBeamElementsAlongLine
 27
 28import numpy as np
 29
 30testIsActive = exu.sys.get('testIsActive', False)
 31
 32OV = exu.OutputVariableType
 33simulationSettings = exu.SimulationSettings()
 34simulationSettings.timeIntegration.newton.useModifiedNewton = False #Just for the test; modified Newton is usually faster
 35testResult = 0
 36
 37def ContactCoordinate(activeConnector):
 38    """a 1D mass at 0.1 above a stop, pushed down by 10 N; returns the system, the node and the contact"""
 39    SC = exu.SystemContainer()
 40    mbs = SC.AddSystem()
 41    nGround = mbs.AddNode(NodePointGround())
 42    node = mbs.AddNode(Node1D(referenceCoordinates=[0], initialCoordinates=[0.1]))
 43    mbs.AddObject(ObjectMass1D(nodeNumber=node, mass=1))
 44    mCoordinate = mbs.AddMarker(MarkerNodeCoordinate(nodeNumber=node, coordinate=0))
 45    mGround = mbs.AddMarker(MarkerNodeCoordinate(nodeNumber=nGround, coordinate=0))
 46    nData = mbs.AddNode(NodeGenericData(numberOfDataCoordinates=1, initialCoordinates=[0.1]))
 47    oContact = mbs.AddObject(ObjectContactCoordinate(markerNumbers=[mGround, mCoordinate], nodeNumber=nData,
 48                                                     contactStiffness=1e4, contactDamping=100,
 49                                                     activeConnector=activeConnector))
 50    mbs.AddLoad(LoadCoordinate(markerNumber=mCoordinate, load=-10))
 51    mbs.Assemble()
 52    mbs.SolveDynamic(simulationSettings)
 53    return (mbs, node, oContact)
 54
 55#(1) contact: at rest on the stop, the gap is the coordinate (stop at 0, no offset): -F/k = -1e-3
 56(mbs, node, oContact) = ContactCoordinate(True)
 57q = mbs.GetNodeOutput(node, OV.Coordinates)  #a float for a node with one coordinate
 58gap = mbs.GetObjectOutput(oContact, OV.Distance)
 59exu.Print('contact: coordinate', q, ', gap', gap, '(both -1e-3 at rest)')
 60testResult += q + gap
 61(mbs, node, oContact) = ContactCoordinate(False)
 62q = mbs.GetNodeOutput(node, OV.Coordinates)  #a float for a node with one coordinate
 63exu.Print('inactive contact: coordinate', q, '(falls freely: 0.1 - F/(2m) t^2 = -4.9)')
 64testResult += q
 65
 66def Planar(addJoint):
 67    """a planar rigid body joined to the ground by addJoint(mbs, mGround, mBody), pushed by a force and a torque"""
 68    SC = exu.SystemContainer()
 69    mbs = SC.AddSystem()
 70    oGround = mbs.AddObject(ObjectGround())
 71    node = mbs.AddNode(NodeRigidBody2D(referenceCoordinates=[0, 0, 0]))
 72    body = mbs.AddObject(ObjectRigidBody2D(nodeNumber=node, mass=2, inertia=0.1))
 73    mGround = mbs.AddMarker(MarkerBodyRigid(bodyNumber=oGround))
 74    mBody = mbs.AddMarker(MarkerBodyRigid(bodyNumber=body))
 75    oJoint = addJoint(mbs, mGround, mBody)
 76    mbs.AddLoad(LoadForceVector(markerNumber=mBody, loadVector=[2, 1, 0]))
 77    mbs.AddLoad(LoadTorqueVector(markerNumber=mBody, loadVector=[0, 0, 0.1]))
 78    mbs.Assemble()
 79    mbs.SolveDynamic(simulationSettings)
 80    return (mbs, node, oJoint)
 81
 82#(2) revolute joint: the joint points stay together
 83(mbs, node, oRevolute) = Planar(lambda mbs, m0, m1: mbs.AddObject(ObjectJointRevolute2D(markerNumbers=[m0, m1])))
 84displacement = mbs.GetObjectOutput(oRevolute, OV.Displacement)
 85exu.Print('revolute joint: displacement', displacement)
 86testResult += np.linalg.norm(displacement)
 87
 88#(3) prismatic joint along x with the rotation constrained: x = F_x/(2m) t^2 = 0.5 at t = 1
 89(mbs, node, oJoint) = Planar(lambda mbs, m0, m1: mbs.AddObject(ObjectJointPrismatic2D(
 90    markerNumbers=[m0, m1], axisMarker0=[2, 0, 0], normalMarker1=[0, 1, 0])))
 91distance = mbs.GetObjectOutput(oJoint, OV.Distance)
 92rotation = mbs.GetObjectOutput(oJoint, OV.Rotation)
 93exu.Print('prismatic joint, rotation constrained: distance', distance, '(0.5), rotation', rotation)
 94testResult += distance + rotation
 95
 96#(4) prismatic joint with the rotation free: the output Rotation is the rotation of the body
 97(mbs, node, oJoint) = Planar(lambda mbs, m0, m1: mbs.AddObject(ObjectJointPrismatic2D(
 98    markerNumbers=[m0, m1], constrainRotation=False)))
 99q = mbs.GetNodeOutput(node, OV.Coordinates)
100rotation = mbs.GetObjectOutput(oJoint, OV.Rotation)
101exu.Print('prismatic joint, rotation free: rotation', rotation, '= body rotation', q[2])
102testResult += rotation + mbs.GetObjectOutput(oJoint, OV.Distance)
103
104#(5) an output variable the item does not declare raises
105refused = 0
106(mbs, node, oRevolute) = Planar(lambda mbs, m0, m1: mbs.AddObject(ObjectJointRevolute2D(markerNumbers=[m0, m1])))
107try:
108    mbs.GetObjectOutput(oRevolute, OV.Rotation)
109except Exception as error:
110    exu.Print('refused:', str(error).split(' [Python file')[0])
111    refused += 1
112
113SC = exu.SystemContainer()
114mbs = SC.AddSystem()
115oGround = mbs.AddObject(ObjectGround())
116beamInfo = GenerateBeamElementsAlongLine(mbs, positionStart=[0, 0, 0], positionEnd=[1, 0, 0],
117                        numberOfElements=1, beamTemplate=ObjectANCFCable2D(massPerLength=1,
118                        bendingStiffness=1, axialStiffness=100), groundConstraintsStart=[1, 1, 1, 1])
119nodes, elements, *_ = beamInfo['nodes'], beamInfo['elements'], beamInfo['loads']
120mCircle = mbs.AddMarker(MarkerBodyPosition(bodyNumber=oGround, localPosition=[0.5, -1, 0]))
121mShape = mbs.AddMarker(MarkerBodyCable2DShape(bodyNumber=elements[0], numberOfSegments=2))
122nData = mbs.AddNode(NodeGenericData(numberOfDataCoordinates=2, initialCoordinates=[0.1, 0.1]))
123oContactCable = mbs.AddObject(ObjectContactCircleCable2D(markerNumbers=[mCircle, mShape], nodeNumber=nData,
124                              numberOfContactSegments=2, circleRadius=0.1, contactStiffness=1e3))
125mbs.Assemble()
126try:
127    mbs.GetObjectOutput(oContactCable, OV.Distance)
128except Exception as error:
129    exu.Print('refused:', str(error).split(' [Python file')[0])
130    refused += 1
131testResult += refused
132
133exu.Print('solution of connectorOutputVariablesTest=', testResult)
134exu.sys['testResult'] = testResult