Cannot Solve Eigenvalues

[utku]

Hi,
When I run my code, I'm having the error below. How can I solve this issue?
Error:

Code: Select all

FullGenEigenSolver::solve() - the eigenvalue 1 is numerically undetermined or infinite
FullGenEigenSolver::solve() - the eigenvalue 2 is numerically undetermined or infinite
FullGenEigenSolver::solve() - the eigenvalue 3 is numerically undetermined or infinite
FullGenEigenSolver::solve() - the eigenvalue 4 is numerically undetermined or infinite
FullGenEigenSolver::solve() - the eigenvalue 5 is numerically undetermined or infinite
FullGenEigenSolver::solve() - the eigenvalue 6 is numerically undetermined or infinite
FullGenEigenSolver::solve() - the eigenvalue 7 is numerically undetermined or infinite
FullGenEigenSolver::solve() - the eigenvalue 8 is numerically undetermined or infinite
FullGenEigenSolver::solve() - the eigenvalue 9 is numerically undetermined or infinite
FullGenEigenSolver::solve() - the eigenvalue 10 is numerically undetermined or infinite

Code:

Code: Select all

import openseespy.opensees as ops
import numpy as np
import matplotlib.pyplot as plt
import opsvis as opsv
# import vfo.vfo as vfo
import more_itertools as mit

ops.wipe()

################################################################################
# Variables
################################################################################

# Bay numbers
x_bay = 4
y_bay = 4
z_bay = 3  # number of storey

# Bay dimensions (m)
x_dist = 5.
y_dist = 4.5
z_dist = 3.2

# Column dimensions (cm => m)
col_x = 40 / 100
col_y = 40 / 100

# Column weight (kN)
col_m = 9.6

# Beam dimensions (cm => m)
beam_w = 40 / 100
beam_h = 60 / 100

# Beam weight (kN)
beam_m = 12.56

# Young modulus of C30/37
# https://eurocodeapplied.com/design/en1992/concrete-design-properties
E = 32837.

# Poisson's ratio
nu = 0.2

# Shear modulus
G = E / (2*(1 + nu))

################################################################################
# Pre-calculations
################################################################################

# Area calculations
a_c = col_x * col_y
a_b = beam_w * beam_h


def torsional_moment(a, b):
    return (a * b**3) * ((16/3) - (3.36 * (b/a) * (1-(b**4/(12*(a**4))))))


def bending_moment(b, h):
    return (b*(h**3))/12


# Torsional moment of inertia calculations
j_c = torsional_moment(max(col_x, col_y), min(col_x, col_y))  # columns
j_b = torsional_moment(max(beam_h, beam_w), min(beam_h, beam_w))  # beams

# Bending moment of inertia calculations
i_y_c = bending_moment(col_x, col_y)  # column, in y axis
i_z_c = bending_moment(col_y, col_x)  # column, in z axis

i_y_b = bending_moment(beam_h, beam_w)  # beam, in y axis
i_z_b = bending_moment(beam_w, beam_h)  # beam, in z axis

################################################################################
# Model Creation
################################################################################

ops.model('basic', '-ndm', 3, '-ndf', 6)

colTransf = 1
beamTransf = 2
girdTransf = 3

coordTransf = 'Linear'
ops.geomTransf(coordTransf, colTransf, 0., -1., 0.)
ops.geomTransf(coordTransf, beamTransf, 0., -1., 0.)
ops.geomTransf(coordTransf, girdTransf, -1., 0., 0.)

################################################################################
# NODES
################################################################################

# Creating coordinates for nodes
coordinates = []
for z in range(0, z_bay+1):
    for y in range(0, y_bay+1):
        for x in range(0, x_bay+1):
            coordinates.append([x*x_dist, y*y_dist, z*z_dist])

# Creating nodes
node_ids = []
for c in range(len(coordinates)):
    node_id = c+1
    ops.node(node_id, coordinates[c][0], coordinates[c][1], coordinates[c][2])
    node_ids.append(node_id)


################################################################################
# COLUMNS
################################################################################

# Creating columns
elem_id = 1
elements = []
for k in range(len(node_ids)):
    top_node = (node_ids[k] + ((x_bay+1)*(y_bay+1)))
    bottom_node = node_ids[k]
    if top_node in node_ids:
        ops.element('elasticBeamColumn', elem_id, bottom_node, top_node,
                    a_c, E, G, j_c, i_y_c, i_z_c, colTransf, col_m)
        elements.append(elem_id)
        elem_id += 1

# Fixing columns
nodes = node_ids.copy()
pieces = (x_bay+1) * (y_bay+1)
nodes = list(mit.chunked(nodes, pieces))

for i in range(0, len(nodes[0])):
    nodeTag = nodes[0][i]
    ops.fix(nodeTag, 1, 1, 1, 1, 1, 1)

################################################################################
# BEAMS
################################################################################

# creating beams on x axis (beams)
for z in range(1, z_bay+1):
    for y in range(y_bay+1):
        for x in range(x_bay):
            bottom_node = x + y*(x_bay+1) + z*(x_bay+1)*(y_bay+1) + 1
            top_node = bottom_node + 1
            ops.element('elasticBeamColumn', elem_id, bottom_node, top_node,
                        a_b, E, G, j_b, i_y_b, i_z_b, beamTransf, beam_m)
            elem_id += 1

# creating beams on y axis (girders)
for z in range(1, z_bay+1):
    for y in range(y_bay):
        for x in range(x_bay+1):
            bottom_node = x + y*(x_bay+1) + z*(x_bay+1)*(y_bay+1) + 1
            top_node = bottom_node + (x_bay+1)
            ops.element('elasticBeamColumn', elem_id, bottom_node, top_node,
                        a_b, E, G, j_b, i_y_b, i_z_b, girdTransf, beam_m)
            elem_id += 1


################################################################################
# RIGID DIAPHRAGM
################################################################################

# TODO: CREATE A FUNCTION TO FIND CENTER OF RIGIDITY
# Creating extra nodes for rigid diaphragm
for i in range(1, z_bay+1):
    node = node_ids[-1]+1
    rigid_x = (x_dist*x_bay)/2
    rigid_y = (y_dist*y_bay)/2
    ops.node(node, rigid_x, rigid_y, z_dist*i)
    node_ids.append(node)

nodes = node_ids.copy()
pieces = (x_bay+1) * (y_bay+1)
nodes = list(mit.chunked(nodes, pieces))

# fixing & binding the nodes for rigid diaphragm
for i in range(len(nodes[-1])):
    ops.fix(nodes[-1][i], 0, 0, 1, 1, 1, 0)
    ops.rigidDiaphragm(3, nodes[-1][i], *nodes[i+1])


################################################################################
# Analysis
################################################################################

ops.constraints('Transformation')
ops.numberer('RCM')
ops.system('BandGeneral')
ops.test('NormDispIncr', 1.0e-6, 6, 2)
ops.algorithm('Linear')
ops.integrator('LoadControl', 1.0)
ops.analysis('Static')
ops.eigen('-fullGenLapack', 10)
ops.analyze(10)

################################################################################
# Plot
################################################################################

opsv.plot_model()

Have a nice day.

[mhscott]

There is an error in your model definition, e.g., insufficient boundary conditions or floating nodes.

[utku]

Do you mean I should fix the 'ops.model()' part or 'ops.fix()' part?

[mhscott]

I see. Your model has no mass.

[utku]

I have defined mass per element as the last parameter of ops.element() command. Should I add another mass with ops.mass() command?

[mhscott]

I didn't see that. You'll have to figure out what's going on. Start with simple models and build up from there.

[utku]

Thank you for your time. Have a nice day.

[ngtranoise]

Thank you for sharing. It looks so troublesome