Natural frequencies are calculated wrong?

Looks that we might have a problem with either stiffness or mass matrix. Here is the test code which can be used to validate the results:

using JuliaFEM
using JuliaFEM.Preprocess
using FEMBase.Test
using Logging
Logging.configure(level=INFO)
add_elements! = JuliaFEM.add_elements!

# Total length of beam
L = 2.0
# Number of elements in beam
nel = 400

# Create uniform 1d mesh for beam elements, first nodes:
X = Dict{Int64, Vector{Float64}}()
for (j, x) in enumerate(linspace(0, L, nel+1))
    X[j] = [x, 0.0, 0.0]
end
nnodes = length(X)
info("Number of nodes: $nnodes")

beam_elements = [Element(Seg2, [j, j+1]) for j=1:nel]
info("Number of elements: ", length(beam_elements))
update!(beam_elements, "geometry", X)
update!(beam_elements, "youngs modulus", 100.0e6)
update!(beam_elements, "shear modulus", 80.0e6)
update!(beam_elements, "density", 100.0)
update!(beam_elements, "cross-section area", 0.1)
update!(beam_elements, "torsional moment of inertia 1", 1.0e-4)
update!(beam_elements, "torsional moment of inertia 2", 1.0e-4)
update!(beam_elements, "polar moment of inertia", 3.0e-4)
update!(beam_elements, "normal", [0.0, 0.0, -1.0])

bc_elements = [Element(Poi1, [1]), Element(Poi1, [3])]
update!(bc_elements, "geometry", X)
for i=1:3
    update!(bc_elements, "fixed displacement $i", 0.0)
    update!(bc_elements, "fixed rotation $i", 0.0)
end

problem = Problem(Beam, "fixed-fixed beam, 2 elements", 6)
add_elements!(problem, beam_elements)
add_elements!(problem, bc_elements)

step = Analysis(Modal)
add_problems!(step, [problem])
run!(step)

freqs = sqrt.(step.properties.eigvals) / (2*pi)

# Calculate accurate solution
E = 100.0e6
I = 1.0E-4
A = 0.1
rho = 100.0
m = rho*A
L = 2.0
acc = 22.3733 * sqrt(E*I/(m*L^4)) / (2*pi)

info("lowest eigenfrequency: ", first(freqs))
info("accurate solution: ", acc)
difference = abs(acc-first(freqs))/max(acc, first(freqs)) * 100
println("difference [%]: ", difference)
@test isapprox(d, 0.0; atol=0.1)

Now, the first natural frequency for fixed-fixed beam should be, according to [1],

With one element, the difference is only 2 % but looks that when the number of elements is increased, the first natural frequency gets lower and lower with our beam formulation. So maybe there is something that depends on the length of the element (for one element, the length is 2, but for two elements the length is 1 and so on).

References

[1] Talukdar S., Vibration of continuous systems. http://www.logosfoundation.org/instrum_gwr/rodo/stalukdar.pdf

JuliaFEM / FEMBeam.jl

Natural frequencies are calculated wrong? #9

References