import numpy as np from scipy.linalg import eigh # Definição das propriedades do material e seção transversal E = 210e9 A = 0.02 Iy = 10e-5 Iz = 20e-5 Gs = 84e9 Jt = 5e-5 rho = 7800 rg = np.sqrt(Iy / A) # Definição da geometria e conectividade dos nós nnode = 8 nel = 8 nnel = 2 ndof = 6 edof = nnel * ndof sdof = nnode * ndof gcoord = np.array([[0, 0, 0], [0, 0, 4], [4, 0, 4], [4, 0, 0], [0, 5, 0], [0, 5, 4], [4, 5, 4], [4, 5, 0]]) # Carregamento aplicado nos nós F = np.array([0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 15000, 0, 0, 0, 0, 0, 0, 0, 0, 0]) nodes = np.array([[0, 4], [1, 5], [2, 6], [3, 7], [4, 5], [5, 6], [6, 7], [7, 4]]) # Definição das restrições nos nós bcdof = np.array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23]) bcval = np.array([0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0, 0]) # Inicialização das matrizes de rigidez e massa globais KK = np.zeros((sdof, sdof)) MM = np.zeros((sdof, sdof)) # Montagem das matrizes de rigidez e massa for i in range(nel): deltax = gcoord[nodes[i, 1], 0] - gcoord[nodes[i, 0], 0] deltay = gcoord[nodes[i, 1], 1] - gcoord[nodes[i, 0], 1] deltaz = gcoord[nodes[i, 1], 2] - gcoord[nodes[i, 0], 2] L = np.sqrt(deltax**2 + deltay**2 + deltaz**2) alpha = np.arctan2(deltay, deltax) beta = np.arctan2(deltaz, np.sqrt(deltax**2 + deltay**2)) c1 = np.cos(alpha) s1 = np.sin(alpha) c2 = np.cos(beta) s2 = np.sin(beta) # Matriz de rotação T T = np.array([[c1, s1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [-s1, c1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, c2, s2, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, -s2, c2, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, c1, s1, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, -s1, c1, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 0, 0, c2, s2, 0], [0, 0, 0, 0, 0, 0, 0, 0, 0, -s2, c2, 0], [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1]]) ke = np.array([[E * A / L, 0, 0, 0, 0, 0, -E * A / L, 0, 0, 0, 0, 0], [0, 12 * E * Iz / L**3, 0, 0, 0, 6 * E * Iz / L**2, 0, -12 * E * Iz / L**3, 0, 0, 0, 6 * E * Iz / L**2], [0, 0, 12 * E * Iy / L**3, 0, -6 * E * Iy / L**2, 0, 0, 0, 6 * E * Iy / L**2, 0, -12 * E * Iy / L**3, 0], [0, 0, 0, Gs * Jt / L, 0, 0, 0, 0, 0, -Gs * Jt / L, 0, 0], [0, 0, -6 * E * Iy / L**2, 0, 4 * E * Iy / L, 0, 0, 0, -6 * E * Iy / L**2, 0, 2 * E * Iy / L, 0], [0, 6 * E * Iz / L**2, 0, 0, 0, 4 * E * Iz / L, 0, -6 * E * Iz / L**2, 0, 0, 0, 2 * E * Iz / L], [-E * A / L, 0, 0, 0, 0, 0, E * A / L, 0, 0, 0, 0, 0], [0, -12 * E * Iz / L**3, 0, 0, 0, -6 * E * Iz / L**2, 0, 12 * E * Iz / L**3, 0, 0, 0, -6 * E * Iz / L**2], [0, 0, 6 * E * Iy / L**2, 0, -6 * E * Iy / L**2, 0, 0, 0, 4 * E * Iy / L, 0, 6 * E * Iy / L**2, 0], [0, 0, 0, -Gs * Jt / L, 0, 0, 0, 0, 0, Gs * Jt / L, 0, 0], [0, 0, -12 * E * Iy / L**3, 0, 2 * E * Iy / L, 0, 0, 0, 6 * E * Iy / L**2, 0, 4 * E * Iy / L, 0], [0, 6 * E * Iz / L**2, 0, 0, 0, 2 * E * Iz / L, 0, -6 * E * Iz / L**2, 0, 0, 0, 4 * E * Iz / L]]) me = (rho * A * L) * np.array([[140 / 420, 0, 0, 0, 0, 0, 70 / 420, 0, 0, 0, 0, 0], [0, (13 / 35) + (6 * rg**2) / (5 * L**2), 0, 0, 0, (11 * L / 210) + (rg**2) / (10 * L), 0, (9 / 70) - (6 * rg**2) / (5 * L**2), 0, 0, (-13 * L / 420) + (rg**2) / (10 * L), 0], [0, 0, (13 / 35) + (6 * rg**2) / (5 * L**2), 0, (-11 * L / 210) - (rg**2) / (10 * L), 0, 0, 0, (9 / 70) - (6 * rg**2) / (5 * L**2), 0, 0, (-13 * L / 420) + (rg**2) / (10 * L)], [0, 0, 0, 140 / 420, 0, 0, 0, 0, 0, 70 / 420, 0, 0], [0, 0, -11 * L / 210 - (rg**2) / (10 * L), 0, (13 / 35) + (6 * rg**2) / (5 * L**2), 0, 0, 0, -13 * L / 420 + (rg**2) / (10 * L), 0, (-9 / 70) - (6 * rg**2) / (5 * L**2), 0], [0, 11 * L / 210 + (rg**2) / (10 * L), 0, 0, 0, (13 / 35) + (6 * rg**2) / (5 * L**2), 0, (-9 / 70) - (6 * rg**2) / (5 * L**2), 0, 0, (-13 * L / 420) + (rg**2) / (10 * L), 0], [70 / 420, 0, 0, 0, 0, 0, 140 / 420, 0, 0, 0, 0, 0], [0, -13 * L / 420 + (rg**2) / (10 * L), 0, 0, 0, (-9 / 70) - (6 * rg**2) / (5 * L**2), 0, (13 / 35) + (6 * rg**2) / (5 * L**2), 0, 0, (11 * L / 210) + (rg**2) / (10 * L), 0], [0, 0, 9 / 70 - (6 * rg**2) / (5 * L**2), 0, -13 * L / 420 + (rg**2) / (10 * L), 0, 0, 0, (13 / 35) + (6 * rg**2) / (5 * L**2), 0, 0, (11 * L / 210) + (rg**2) / (10 * L)], [0, 0, 0, 70 / 420, 0, 0, 0, 0, 0, 140 / 420, 0, 0], [0, -13 * L / 420 + (rg**2) / (10 * L), 0, 0, (-9 / 70) - (6 * rg**2) / (5 * L**2), 0, 0, (11 / 35) + (6 * rg**2) / (5 * L**2), 0, 0, (11 * L / 210) + (rg**2) / (10 * L), 0], [0, 0, -13 * L / 420 + (rg**2) / (10 * L), 0, 0, (-13 * L / 420) + (rg**2) / (10 * L), 0, 0, (-11 * L / 210) - (rg**2) / (10 * L), 0, 0, (L**2 / 105) + (2 * rg**2) / (15)]]) Ke = np.dot(np.dot(T.T, ke), T) Me = np.dot(np.dot(T.T, me), T) idx = np.zeros((1, edof), dtype=int) for j in range(nnel): idx[0, j * ndof:(j + 1) * ndof] = np.arange(ndof) + nodes[i, j] * ndof for m in range(edof): for n in range(edof): KK[idx[0, m], idx[0, n]] += Ke[m, n] MM[idx[0, m], idx[0, n]] += Me[m, n] # Redução das matrizes devido às condições de contorno KK_cc = np.delete(np.delete(KK, bcdof, axis=0), bcdof, axis=1) MM_cc = np.delete(np.delete(MM, bcdof, axis=0), bcdof, axis=1) F_cc = np.delete(F, bcdof) print("Dimensão de KK_cc:", KK_cc.shape) print("Dimensão de F_cc:", F_cc.shape) # Resolução do sistema de equações lineares u = np.linalg.solve(KK_cc, F_cc) # Recuperação dos deslocamentos totais U = np.zeros((1, sdof)) idx_cc = 0 idx_gdl = 0 for i in range(sdof): if idx_cc < len(bcdof) and i == bcdof[idx_cc]: U[0, i] = bcval[idx_cc] idx_cc += 1 else: U[0, i] = u[idx_gdl] idx_gdl += 1 Ut = np.transpose(U) # Cálculo das reações de apoio reacao = np.dot(KK, Ut) # Determinação das frequências naturais e modos de vibração W, V = eigh(KK_cc, MM_cc) f = (np.real(W))**0.5 / (2*np.pi) # Impressão dos resultados print("Deslocamentos totais:") print(U) print("\nReações de apoio:") print(reacao) print("\nFrequências naturais:") for i in range(20): print('{:.3f} Hz'.format(f[i])) print("\nModos de vibração correspondentes:") print(np.real(V))