I have a few questions on the strategy for my problem setup.

My problem is made of a group containing 5 key components, but 4 are instances of the same component (cable), which calculates the tension (its 3 components in 3D space) based on the deflections uvw (array(3)) at one end of the cable. The last component (balance) contains a nonlinear system of equations, which takes external loads as parameters, and the cable tensions as inputs; it outputs the deflections at a common point uvw_tip (array(3)), theta_tip, and the reactions XYZ0 (array(3)). A total of 7 output scalars.

The outputs uvw_tip and theta_tip are processed through an auxiliary component to calculate the uvw to input back into the cable components.

So I believe my balance component is an ImplicitComponent, yet the theta_tip does not even show up directly in the residual equations, because its effect is buried in the other inputs from the cable components.

There is also the coupling issue of connecting this to cable, and I am not sure of the best approach here, and not sure what I have setup so far is appropriate. If you could recommend a direction, I would appreciate it.

A related question is: in the ImplicitComponent I have 7 residual equations, but I am not interested in all the state variables, and only in the outputs mentioned above.
How do I handle the fact that the uvw_tip and XYZ0 are arrays in the residuals[xxx] below? Thank you in advance.

residuals['theta_legtip']  = L_L0**2 - ((L_L0  + uvw_legtip[0])**2 + uvw_legtip[1]**2 + uvw_legtip[2]**2)  #length of leg is invariant   ::Eq#1
residuals['XYZ0'][0] = Fxyz[0] + XYZ0[0] + Nxyz_UC1[0] + Nxyz_UC2[0] + Nxyz_LC1[0] + Nxyz_LC2[0] #Balance along X                ::Eq#2  
residuals['XYZ0'][1] = Fxyz[1] + XYZ0[1] + Nxyz_UC1[1] + Nxyz_UC2[1] + Nxyz_LC1[1] + Nxyz_LC2[1] #Balance along Y                ::Eq#3  
residuals['XYZ0'][2] = Fxyz[2] + XYZ0[2] + Nxyz_UC1[2] + Nxyz_UC2[2] + Nxyz_LC1[2] + Nxyz_LC2[2] + M_L*gravity #Balance along Z  ::Eq#4  
residuals['uvw_legtip'][0] = Mxyz[0] - L_L[2]*Fxyz[1] + (Fxyz[2] - M_L * gravity/2.) * L_L[1]          + 
                                                              UC1_Bp[1]*Nxyz_UC1[2] - UC1_Bp[2]*Nxyz_UC1[1] + 
                                                              UC2_Bp[1]*Nxyz_UC2[2] - UC2_Bp[2]*Nxyz_UC2[1] + 
                                                              LC1_Bp[1]*Nxyz_LC1[2] - LC1_Bp[2]*Nxyz_LC1[1] + 
                                                              LC2_Bp[1]*Nxyz_LC2[2] - LC2_Bp[2]*Nxyz_LC2[1]   #Eq.#5
residuals['uvw_tip'][1] = Mxyz[1] + L_L[1]*Fxyz[0] + (-Fxyz[2] + M_L * gravity/2.) * L_L[0]          + 
                                                              - UC1_Bp[0]*Nxyz_UC1[2] + UC1_Bp[2]*Nxyz_UC1[0]  
                                                              - UC2_Bp[0]*Nxyz_UC2[2] + UC2_Bp[2]*Nxyz_UC2[0]  
                                                              - LC1_Bp[0]*Nxyz_LC1[2] + LC1_Bp[2]*Nxyz_LC1[0]  
                                                              - LC2_Bp[0]*Nxyz_LC2[2] + LC2_Bp[2]*Nxyz_LC2[0]  #Eq.#6
        
residuals['uvw_tip'][2] = Mxyz[2] + L_L[0]*Fxyz[1] - L_L[1]*Fxyz[0]                                  + 
                                                              - UC1_Bp[1]*Nxyz_UC1[0] + UC1_Bp[0]*Nxyz_UC1[1]  
                                                              - UC2_Bp[1]*Nxyz_UC2[0] + UC2_Bp[0]*Nxyz_UC2[1]  
                                                              - LC1_Bp[1]*Nxyz_LC1[0] + LC1_Bp[0]*Nxyz_LC1[1]  
                                                              - LC2_Bp[1]*Nxyz_LC2[0] + LC2_Bp[0]*Nxyz_LC2[1]  #Eq.#7