@@ -3895,8 +3895,8 @@ void MultibodyPlant<T>::CalcReactionForces(
3895
3895
// for the discrete solvers we have today, though it might change for future
3896
3896
// solvers that prefer an implicit evaluation of these terms.
3897
3897
// TODO(amcastro-tri): Consider having a
3898
- // DiscreteUpdateManager::EvalReactionForces() to ensure the manager performs
3899
- // this computation consistently with its discrete update.
3898
+ // DiscreteUpdateManager::EvalReactionForces() to ensure the manager performs
3899
+ // this computation consistently with its discrete update.
3900
3900
const VectorX<T>& vdot = this ->EvalForwardDynamics (context).get_vdot ();
3901
3901
std::vector<SpatialAcceleration<T>> A_WB_vector (num_bodies ());
3902
3902
std::vector<SpatialForce<T>> F_BMo_W_vector (num_bodies ());
@@ -3908,9 +3908,9 @@ void MultibodyPlant<T>::CalcReactionForces(
3908
3908
// Since vdot is the result of Fapplied and tau_applied we expect the result
3909
3909
// from inverse dynamics to be zero.
3910
3910
// TODO(amcastro-tri): find a better estimation for this bound. For instance,
3911
- // we can make an estimation based on the trace of the mass matrix (Jain
3912
- // 2011, Eq. 4.21). For now we only ASSERT though with a better estimation we
3913
- // could promote this to a DEMAND.
3911
+ // we can make an estimation based on the trace of the mass matrix (Jain 2011,
3912
+ // Eq. 4.21). For now we only ASSERT though with a better estimation we could
3913
+ // promote this to a DEMAND.
3914
3914
// TODO(amcastro-tri) Uncomment this line once issue #12473 is resolved.
3915
3915
// DRAKE_ASSERT(tau_id.norm() <
3916
3916
// 100 * num_velocities() *
@@ -3924,49 +3924,61 @@ void MultibodyPlant<T>::CalcReactionForces(
3924
3924
internal_tree ().get_joint_mobilizer (joint_index);
3925
3925
const internal::Mobilizer<T>& mobilizer =
3926
3926
internal_tree ().get_mobilizer (mobilizer_index);
3927
-
3928
- // Reversed means the joint's parent(child) body is the outboard(inboard)
3929
- // body for the mobilizer.
3930
- const bool is_reversed = mobilizer.mobod ().is_reversed ();
3927
+ const internal::MobodIndex mobod_index = mobilizer.mobod ().index ();
3931
3928
3932
3929
// F_BMo_W is the mobilizer reaction force on mobilized body B at the origin
3933
3930
// Mo of the mobilizer's outboard frame M, expressed in the world frame W.
3934
- const SpatialForce<T>& F_BMo_W = F_BMo_W_vector[mobilizer_index];
3935
-
3936
- // But the quantity of interest, F_CJc_Jc, is the joint's reaction force on
3937
- // the joint's child body C at the joint's child frame Jc, expressed in Jc.
3938
- SpatialForce<T>& F_CJc_Jc = output->at (joint.ordinal ());
3931
+ const SpatialForce<T>& F_BMo_W = F_BMo_W_vector[mobod_index];
3939
3932
3940
3933
// Frames of interest:
3934
+ const Frame<T>& frame_Jp = joint.frame_on_parent ();
3941
3935
const Frame<T>& frame_Jc = joint.frame_on_child ();
3942
- const Frame<T>& frame_M = mobilizer.outboard_frame ();
3936
+ const FrameIndex F_index = mobilizer.inboard_frame ().index ();
3937
+ const FrameIndex M_index = mobilizer.outboard_frame ().index ();
3938
+ const FrameIndex Jp_index = frame_Jp.index ();
3939
+ const FrameIndex Jc_index = frame_Jc.index ();
3940
+
3941
+ // In Drake we must have either:
3942
+ // - Jp == F and Jc == M (typical case)
3943
+ // - Jp == M and Jc == F (mobilizer is reversed from joint)
3944
+ DRAKE_DEMAND ((Jp_index == F_index && Jc_index == M_index) ||
3945
+ (Jp_index == M_index && Jc_index == F_index));
3946
+
3947
+ // Mobilizer is reversed if the joint's parent frame Jp is the mobilizer's
3948
+ // outboard frame M.
3949
+ const bool is_reversed = (Jp_index == M_index);
3943
3950
3944
3951
// We'll need this in both cases below since we're required to report
3945
3952
// the reaction force expressed in the joint's child frame Jc.
3946
3953
const RotationMatrix<T> R_JcW =
3947
3954
frame_Jc.CalcRotationMatrixInWorld (context).inverse ();
3948
3955
3949
- if (&frame_M == &frame_Jc) {
3950
- // This is the easy case. Just need to re-express.
3951
- F_CJc_Jc = R_JcW * F_BMo_W;
3952
- continue ;
3953
- }
3956
+ // The quantity of interest, F_CJc_Jc, is the joint's reaction force on the
3957
+ // joint's child body C at the joint's child frame Jc, expressed in Jc.
3958
+ SpatialForce<T>& F_CJc_Jc = output->at (joint.ordinal ());
3959
+ if (!is_reversed) {
3960
+ F_CJc_Jc = R_JcW * F_BMo_W; // The typical case: Mo==Jc and B==C.
3961
+ } else {
3962
+ // For this reversed case, F_BMo_W is the reaction on the joint's _parent_
3963
+ // body at Jp, expressed in W.
3964
+ const SpatialForce<T>& F_PJp_W = F_BMo_W; // Reversed: Mo==Jp and B==P.
3954
3965
3955
- // If the mobilizer is reversed, Newton's 3ʳᵈ law (action/reaction) (and
3956
- // knowing Drake's joints are massless) says the force on the child at M
3957
- // is equal and opposite to the force on the parent at M .
3958
- const SpatialForce<T> F_CMo_W = is_reversed ? -F_BMo_W : F_BMo_W ;
3959
- const SpatialForce<T> F_CMo_Jc = R_JcW * F_CMo_W ; // Reexpress in Jc.
3966
+ // Newton's 3ʳᵈ law (action/reaction) (and knowing Drake's joints are
3967
+ // massless) says the force on the child _at Jp_ is equal and opposite to
3968
+ // the force on the parent at Jp .
3969
+ const SpatialForce<T> F_CJp_W = -F_PJp_W ;
3970
+ const SpatialForce<T> F_CJp_Jc = R_JcW * F_CJp_W ; // Reexpress in Jc.
3960
3971
3961
- // However, the reaction force we want to report on the child is at Jc,
3962
- // not M . We need to shift the application point from Mo to Jco .
3972
+ // However, the reaction force we want to report on the child is at Jc,
3973
+ // not Jp . We need to shift the application point from Jp to Jc .
3963
3974
3964
- // Find the shift vector p_MoJco_Jc (= -p_JcoMo_Jc ).
3965
- const RigidTransform<T> X_JcM = frame_M .CalcPose (context, frame_Jc);
3966
- const Vector3<T> p_MoJco_Jc = -X_JcM .translation ();
3975
+ // Find the shift vector p_JpJc_Jc (= -p_JcJp_Jc ).
3976
+ const RigidTransform<T> X_JcJp = frame_Jp .CalcPose (context, frame_Jc);
3977
+ const Vector3<T> p_JpJc_Jc = -X_JcJp .translation ();
3967
3978
3968
- // Perform the M->Jc shift.
3969
- F_CJc_Jc = F_CMo_Jc.Shift (p_MoJco_Jc);
3979
+ // Perform the Jp->Jc shift.
3980
+ F_CJc_Jc = F_CJp_Jc.Shift (p_JpJc_Jc);
3981
+ }
3970
3982
}
3971
3983
}
3972
3984
0 commit comments