Compliance Not Considering Point Loads

[JohnDN90]

It appears that compliance_output does not consider point loads when it is being evaluated here.

I tested it on 2D representation of a cantilever beam, bending in-plane, both with 1 element and 4 elements in the mesh. For both cases, the vec vector remains entirely zero after e.side_external_residual and e.volume_external_residual have been called.

I only started diving into the "output" type source code today, but to me it looks like we'll need a way to determining which node(s) that point load is acting and which element(s) belong to that node. Thenm apply it in a way that the point load doesn't get redundantly included if it acts on a node that is shared by multiple elements.

Do you have any suggestions/thoughts on whats the best/easiest way to include point loads in the compliance output?

[manavbhatia]

Point loads in structural analysis always make for interesting cases.

Output base class in MAST is currently based on the notion of element operations. I don't think putting a hack inside the output class with the nodes would be a good idea. This should be done at a higher level in the assembly routines.

I think there are a couple different options to address this:

1. libMesh defines a NodeElem that are wrappers around nodes so that node-level operations can be fit within the same context as elements. This way the existing element iterators used in the assembly routines also accommodate for nodal quantities. For example, you could define NodeElems for each node on which point loads are applied and then associated your SPC and point load data on node elements, instead of defining them on nodes. Then, the assembly routines can be setup to also include these elements in their global assembly routines. This would require some modest amount of code modification outside of Output class.

2. Currently the output assembly done here only iterates over the elements and does not include nodal quantities. Likewise for the following routines that compute the RHS for adjoint solution and sensitivity analysis. For each of these cases we could add a for loop at the end that operates on the nodal quantities, as necessary. We have done something similar for the nodal loads in the nonlinear analysis here.

[JohnDN90]

Thanks for the input. I was messing around with putting a hack inside the the output class and while it worked when running with 1 processor, I ran into issues with multiple processors.

I think option 2 that you presented might be better than option 1 for a couple reasons. First, I think the concept of a node acting as an element isn't very intuitive and might be a bit confusing for adopters of MAST who have a commercial FEA program background.

Second, Josh and I have discussed adding in zero length elements in the future using the NodeElem in libMesh to model the CBUSH or CELAS elements with coincident nodes in Nastran. (We currently have zero length bushings supported in one our branches, but are using a hack to perturb 1 node by a very small amount so that the nodes are no longer coincident.) If we do proceed with that, I think having NodeElems that represent elements and NodeElems that are just nodes that need to act like elements would be a bit confusing as well.

[JohnDN90]

A couple more quick questions on compliance_output.

It looks like here you're making the assumption that the side_external_residual and the volume_external_residual do not depend on the displacements. This would mean this formulation for compliance isn't valid for follower loads. Is my understanding correct?

Also, if I recall correctly, you are considering the thermoelastic forces part of the internal residual, instead of the external residual. Is that correct? If this is the case, it looks like thermoelastic loads are going to be ignored in the calculation of compliance.

[manavbhatia]

Second, Josh and I have discussed adding in zero length elements in the future using the NodeElem in libMesh to model the CBUSH or CELAS elements with coincident nodes in Nastran. (We currently have zero length bushings supported in one our branches, but are using a hack to perturb 1 node by a very small amount so that the nodes are no longer coincident.) If we do proceed with that, I think having NodeElems that represent elements and NodeElems that are just nodes that need to act like elements would be a bit confusing as well.

What is the purpose of CBUSH and CELAS? If it truly is a single node element then I would not recommended adding two separate nodes. If the element Jacobian is not going to be computed then it may be best to specify the same node for both ends of the elements.

I think it might be worth looking into the NodeElem for this perspective. I recognize that this may be confusing to folks more ingrained in Nastran terminology, but I suspect that the creation of NodeElems would all happen behind the scenes and they would not have to directly engage with them.

[manavbhatia]

It looks like here you're making the assumption that the side_external_residual and the volume_external_residual do not depend on the displacements. This would mean this formulation for compliance isn't valid for follower loads. Is my understanding correct?

No.

The element is provided the current solution (displacement) vector before these calls are made. The element routines that compute the residual vector contributions from side and volume loads will be responsible for using the displacement, if that is required for accounting for the nonlinearities.

I should note that we currently do not have an implementation of the follower force, although we can update our residual computation routines (from surface pressure/traction) to account for this.

Also, if I recall correctly, you are considering the thermoelastic forces part of the internal residual, instead of the external residual. Is that correct? If this is the case, it looks like thermoelastic loads are going to be ignored in the calculation of compliance.

The thermoelastic forces are computed in the thermal_residual routine, which is called by volume_external_residual .

[JohnDN90]

Thanks for clarifying on the thermoelastic loads.

The element is provided the current solution (displacement) vector before these calls are made. The element routines that compute the residual vector contributions from side and volume loads will be responsible for using the displacement, if that is required for accounting for the nonlinearities.

I'm still not quite seeing it. Let C be compliance, p be a general parameter, f be the external forces, u be the displacements. d denotes the total derivative while dp denotes the partial derivative. It looks to me like you are calculating the sensitivity as:

dC_dp = dpf_dpp.dot(u) + f.dot(du_dp)

In the general sense, the equation should be:
dC_dp = dpC_dpp + dpC_dpu.dot(du_dp)
where
dpC_dpp = dpf_dpp.dot(u) and dpC_dpu = dpf_dpu.dot(u) + f

I think the dpf_dpu.dot(u) term is missing. If I'm just glancing over it, could you point out in the code where that part is calculated?

[JohnDN90]

What is the purpose of CBUSH and CELAS? If it truly is a single node element then I would not recommended adding two separate nodes. If the element Jacobian is not going to be computed then it may be best to specify the same node for both ends of the elements.

CBUSH and CELAS are essentially spring elements. I'm relatively new to them (in Nastran) myself, but my understanding is that CBUSH is a newer version of CELAS. CELAS doesn't take the orientation between the two elements into account and can introduce additional constraints into the model (I haven't experience this myself, I've just found it mentioned online in forums and a PowerPoint). CBUSH does take the orientation into account. Additionally, unlike the traditional spring element in FEA, the CBUSH element takes length into account. It's stiffness matrix is quite similar to that of a frame element including torsion, extension, and Euler-Bernoulli bending, with the geometric/material stiffnesses replaced by spring constants.

If the element has only 1 node, it acts as a ground spring. With two nodes (non-coincident) it acts as a spring (well, more like a frame as described above). The element can also have two nodes that are coincident. That can be used to represent a variety of things, but I find it easiest to envision it as an adhesive holding two parts together at that node with a stiffness in each direction. The main issue that came out of this was that there wasn't already a way in MAST to specify both the elemental x-axis and y-axis directions, which is necessary since you can't geometrically obtain the element x-axis in that case. Also, if I recall correctly, I think there may have been an issue with libMesh removing one of the points since they were in the same geometric location, but I may be remembering wrong or it could have been a bug I introduced.

So, in order to be able to geometrically obtain the x-axis. I read in the true x-axis from the coordinate system assigned to the CBUSH in the BDF, and then perturb one of the nodes in that direction by a very small amount (~ 1e-07).

So yes, it essentially a hack to support CBUSH with coincident nodes but the error it introduces seems to be negligible.

[manavbhatia]

In the general sense, the equation should be:
dC_dp = dpC_dpp + dpC_dpu.dot(du_dp)
where
dpC_dpp = dpf_dpp.dot(u) and dpC_dpu = dpf_dpu.dot(u) + f

I think the dpf_dpu.dot(u) term is missing. If I'm just glancing over it, could you point out in the code where that part is calculated?

You are right. This needs to be added. However, I think the missing term should be dpf_dpu.transpose().dot(u).

This will need correction at a couple of places in compliance_output.cpp. One in the output_derivative_for_elem() and the other should be in evaluate_sensitivity(). We are currently passing a dummy matrix for the element jacobian. We can instead ask for the Jacobian to be computed and then add the dot product of this Jacobian-transposed with the element solution to the sensitivity/derivative.

[manavbhatia]

CBUSH and CELAS are essentially spring elements. I'm relatively new to them (in Nastran) myself, but my understanding is that CBUSH is a newer version of CELAS. CELAS doesn't take the orientation between the two elements into account and can introduce additional constraints into the model (I haven't experience this myself, I've just found it mentioned online in forums and a PowerPoint). CBUSH does take the orientation into account. Additionally, unlike the traditional spring element in FEA, the CBUSH element takes length into account. It's stiffness matrix is quite similar to that of a frame element including torsion, extension, and Euler-Bernoulli bending, with the geometric/material stiffnesses replaced by spring constants.

Why is length needed if a spring constant is used? Do you specify independent spring constants in the three translation and three rotation axes?

The main issue that came out of this was that there wasn't already a way in MAST to specify both the elemental x-axis and y-axis directions, which is necessary since you can't geometrically obtain the element x-axis in that case.

Yes, the class hierarchy was designed with the intent to use beams/plates/solid elements. So, the use of local_elem becomes relevant, which requires a way to compute the element x- and y- axes.

Also, if I recall correctly, I think there may have been an issue with libMesh removing one of the points since they were in the same geometric location, but I may be remembering wrong or it could have been a bug I introduced.

If you compile libMesh with --enable-unique-id then it will allow coincident nodes to be used. Otherwise it removes coincident nodes.

So, in order to be able to geometrically obtain the x-axis. I read in the true x-axis from the coordinate system assigned to the CBUSH in the BDF, and then perturb one of the nodes in that direction by a very small amount (~ 1e-07).

So, you don't need the element length?

[JohnDN90]

You are right. This needs to be added. However, I think the missing term should be dpf_dpu.transpose().dot(u)

You're correct.

And I've been working on adding a few new output classes, so I'm starting to become familiar with how you have it structured. I should be able to make the necessary changes to include that term fairly easily.

Why is length needed if a spring constant is used?

The CBUSH supports assigning a spring constant for all six degrees of freedom. In older Nastran-based programs, (i.e. MYSTRAN) the CBUSH acts purely as a spring element, where there is no coupling between bending and displacement like there is in beam elements (i.e. y-displacement and bending about the z-axis) and the length is not required. However, in more modern versions, MSC Nastran and Optistruct (not sure about other versions of Nastran), they incorporate this bending-displacement coupling which is dependent on length and the relevant DOF spring constant. Abaqus does the same thing with their element that's equivalent to a CBUSH, the formulation is a little different, but it still depends on length.

Do you specify independent spring constants in the three translation and three rotation axes?

Yes.

So, you don't need the element length?

When nodes are not coincident, you do need length. For CBUSH with coincident nodes, the length is zero and the terms coupling the bending and displacement disappear. In the 'hack' I've implemented, the coupling terms are not eliminated because the length is nonzero. However, the length is very small, so while the coupling terms do still exist, they are small relative to the rest of the terms in the stiffness matrix and have a negligible affect on the analysis.

[JohnDN90]

If you compile libMesh with --enable-unique-id then it will allow coincident nodes to be used. Otherwise it removes coincident nodes.

Just commenting about this for future reference for myself. I compiled libMesh with --disable-unique-id, so that was the issue with my coincident node being removed. I should compile with --enable-unique-id in the future.