Hi everyone,
I am researching a composite beam modeling approach by comparing two cantilever setups:
-
A baseline beam modeled entirely with C3D20 (20-node hexahedral) elements.
-
A composite beam where S8 shell elements are “jacketed” onto the top and bottom surfaces of a C3D20 solid core.
The Setup:
-
Loading: Surface traction on the free end (right face).
-
Boundary Condition: The left face of the solid core is fixed using *BOUNDARY.
-
Interaction: I used *TIE to bind the shell-solid interfaces, with the solid face as the Master and the shell surface as the Slave.
The Issue:
In the initial run, the shells “peeled” away from the solid core specifically at the fixed end. Even though the Master (solid) nodes were fixed, the Slave (shell) nodes at that specific boundary seemed to ignore the constraint.
Attempted Fix:
I applied an additional *BOUNDARY to the fixed-end edges of the shells. This stopped the peeling, but the solver issued warnings indicating that tied constraints were not generated for those nodes because they were already subject to an SPC. However-
The “Paradox” in the Results:
-
Static Analysis: Despite the deactivation of the tie at the fixed boundary, the displacement results for both the pure solid and composite beams match almost perfectly.
-
Modal/Matrix Analysis: When I extract the final global stiffness matrix (I wrote a custom code to export the stiffness vals at the end in mafillsm.f) and solved for the lowest positive eigenvalues, the results revealed that the composite beams were much softer than the pure solid one. The eigenvalues were consistently ~45% lower than those of the pure solid cantilever beam.
Later, I was able to overcome the over constrained nodes issue by swapping the Master and slave nodes. But the low eigenvalue issue still remains.
I now suspect that the stiffness matrix I’m exporting is the root cause of the problem. For the following reasons-
· The issue exists only in the stiffness matrix analysis. The simulation results match well.
· The size of the stiffness matrix for the composite beam scenario is much higher than the pure solid beam scenario. I think if the MPCs were applied correctly, it should be a lot lower as the solid and shell DOFs (the displacements) would be tied by a single equation per node. Unless extra penalty equations are introduced. [I’m looking into this]
It would be great if someone could enlighten me on this as I’m relatively new to numerical code dev.
PS- This is the block of code that exports the K matrix. It is at the bottom of mafillsm.f currently.
c — Exporting the FINAL assembled matrix —
open(unit=99,file='K_final.mtx',status='unknown')
do i=1,neq(1)
write(99,'(2I10,1X,E20.12)') i,i,ad(i)
do j=jq(i),jq(i+1)-1
write(99,'(2I10,1X,E20.12)') i,irow(j),au(j)
end do
end do
close(99)
return
end