Finite element method Beginning 2 (Preparation of one-dimensional overall stiffness matrix)

About this article

Previous article was up to the program for creating the element stiffness matrix, so this time I will explain the program for creating the overall stiffness matrix from the element stiffness matrix. .. Similar to the previous time, it is based on The first finite element method by Java published by Corona Publishing.

One-dimensional finite element method (fabrication of total stiffness matrix)

Combine element stiffness matrices to create an overall stiffness matrix

Last time, I made an element stiffness matrix and an element stiffness equation for one element of an object. We will formulate this for each element in the object. If an element i has nodes j, k, then the element stiffness matrix and the element stiffness equations are:


  {\bf K}^{(i)} \left(
    \begin{array}{c}
      u_j \\
      u_k
    \end{array}
  \right)
  = 
  \left(
    \begin{array}{c}
      K_{11}^{(i)}\ K_{12}^{(i)} \\
      K_{21}^{(i)}\ K_{22}^{(i)}
    \end{array}
  \right)
   \left(
    \begin{array}{c}
      u_j \\
      u_k
    \end{array}
  \right)
  =
  \left(
    \begin{array}{c}
      f^{(i)}_j \\
      f^{(i)}_k
    \end{array}
  \right)
  (1)

Can be written as. For example, if there are three elements and the i-th element has two nodes (i, i + 1), the element stiffness equation is as follows.


  \left(
    \begin{array}{c}
      K_{11}^{(1)}\ K_{12}^{(1)} \\
      K_{21}^{(1)}\ K_{22}^{(1)}
    \end{array}
  \right)
   \left(
    \begin{array}{c}
      u_1 \\
      u_2
    \end{array}
  \right)
  =
  \left(
    \begin{array}{c}
      f^{(1)}_1 \\
      f^{(1)}_2
    \end{array}
  \right)
  (2)


  \left(
    \begin{array}{c}
      K_{11}^{(2)}\ K_{12}^{(2)} \\
      K_{21}^{(2)}\ K_{22}^{(2)}
    \end{array}
  \right)
   \left(
    \begin{array}{c}
      u_2 \\
      u_3
    \end{array}
  \right)
  =
  \left(
    \begin{array}{c}
      f^{(2)}_2 \\
      f^{(2)}_3
    \end{array}
  \right)
  (3)


  \left(
    \begin{array}{c}
      K_{11}^{(3)}\ K_{12}^{(3)} \\
      K_{21}^{(3)}\ K_{22}^{(3)}
    \end{array}
  \right)
   \left(
    \begin{array}{c}
      u_3 \\
      u_4
    \end{array}
  \right)
  =
  \left(
    \begin{array}{c}
      f^{(3)}_3 \\
      f^{(3)}_4
    \end{array}
  \right)
  (4)

Considering that each element shares each node, these equations can be combined into one,


  \left(
    \begin{array}{cccc}
      K_{11}^{(1)}& K_{12}^{(1)} & 0 & 0 \\
      K_{21}^{(1)}& K_{22}^{(1)} + K_{11}^{(2)} & K_{12}^{(2)} & 0 \\
     0 & K_{21}^{(2)} & K_{22}^{(2)} + K_{11}^{(3)} & K_{12}^{(3)} \\
     0 & 0  & K_{21}^{(3)} & K_{22}^{(3)} \\
    \end{array}
  \right)
   \left(
    \begin{array}{c}
      u_1 \\
      u_2 \\
      u_3 \\
      u_4 
    \end{array}
  \right)
  =
  \left(
    \begin{array}{c}
      f^{(1)}_1\\
      f^{(1)}_2+f^{(2)}_2 \\
      f^{(2)}_3+f^{(3)}_3 \\
      f^{(3)}_4
    \end{array}
  \right)
  (5)

It will be. The sum of the element stiffness equations in this way is called the overall stiffness equation, and the coefficient matrix on the left side is called the overall stiffness matrix. The total stress and strain are calculated by solving the overall stiffness equation in consideration of external forces and boundary conditions.

Calculation program for the overall stiffness matrix of one-dimensional elements

When the node coordinates, element area / Young's modulus, force applied to the node, and node displacement constraint condition of each node are input, the program that outputs the corresponding overall stiffness equation is as follows. (The force applied to the nodes and the node displacement constraints of each node are not used this time, but they will be necessary parameters when solving the next overall stiffness equation.) In the program, the Element class created last time is used. I will.


import numpy as np

    #Data entry and calculation preparation
    #list_x :Node coordinates
    #list_area :Area of each element
    #list_young :Young's modulus of each element
    #list_force :Force applied to each node
    #list_nodeCondition :Nodal displacement constraint condition of each node
    def __init__(self,list_x,list_area,list_young,list_force,list_nodeCondition):
        #Setting the number of nodes NumberOfNode, nnode
        self.NumberOfNode = int(list_x.shape[0])
        nnode = self.NumberOfNode
        
        #Setting the number of elements NumberOfElement, nelm
        self.NumberOfElement = nnode-1
        nelm = self.NumberOfElement
        
        #Create lists and arrays for data storage, calculation, and output
        self.x = np.zeros(nnode)                     #Nodal coordinates
        self.force = np.zeros(nnode)                 #Force on the node
        self.nodeCondition = np.zeros((nnode,2))     #Nodal displacement constraint condition
        self.totalStiffness = np.zeros((nnode,nnode))#Overall stiffness matrix
        self.disp = np.zeros(nnode)                  #Nodal displacement
        self.stress = np.zeros(nelm)                 #stress
        
        #List for storing Element classes
        self.elem = np.array([])
        
        #Initial condition storage of node coordinates, force applied to nodes, and node displacement conditions
        self.x = list_x
        self.force = list_force
        self.nodeCondition = list_nodeCondition
        
        #Create Element class and store in elem list
        #Then store the data about the element. Implemented for the number of elements.
        for i in range(nelm):
            self.elem = np.append(self.elem, Element())
            self.elem[i].young=list_young[i]
            self.elem[i].area=list_area[i]
            for j in range(2):
                #Stores the contact numbers ix and coordinates xe that make up each element
                self.elem[i].ix[j]=i+j
                self.elem[i].xe[j]=list_x[i+j]
                
    #Fabrication of total stiffness matrix
    def assembleMatrix(self):
        #Bmatrix / element stiffness matrix creation for each element
        print("\nAssemble Stiffness Matrix : \n")
        nelm = self.NumberOfElement
        for i in range(nelm):
            self.elem[i].makeBmatrix()
            self.elem[i].makeElementStiffness()
        
        #Create an overall stiffness matrix from each element stiffness matrix
        for k in range(nelm):
            print("  Element",k)
            for i in range(2):
                for j in range(2):
                    ii = int(self.elem[k].ix[i])
                    jj = int(self.elem[k].ix[j])
                    self.totalStiffness[ii,jj] += self.elem[k].Kmatrix[i,j] 

The above operation example is as follows. Here, it is assumed that there are four nodes, each coordinate is (50, 150, 250, 350) (mm), the cross-sectional area of each element is 100 mm $ ^ 2 $, and Young's modulus is 200 GPa. Also, although not used this time, list_nodeCondition specifies that the node coordinates are fixed at 50 mm, and list_force specifies that a force of 1000 N can be applied to the node coordinates of 350 mm.

input

#node parameters
list_x = np.array([50,150,250,350])
list_force = np.array([0.0,0.0,0.0,1000.])
list_nodeCondition = np.array([[1,0],[0,0],[0,0],[0,0]])

#Element parameters
list_area = np.array([100,100,100])
list_young = np.array([200000,200000,200000])
 
Fem = Fem1dim(list_x,list_area,list_young,list_force,list_nodeCondition)
Fem.assembleMatrix()
print(Fem.totalStiffness)

output


Assemble Stiffness Matrix : 

  Element 0
  Element 1
  Element 2
[[ 200000. -200000.       0.       0.]
 [-200000.  400000. -200000.       0.]
 [      0. -200000.  400000. -200000.]
 [      0.       0. -200000.  200000.]]

Summary

This time, we explained the formula of the overall stiffness equation from the element stiffness equation and introduced the overall stiffness matrix fabrication program. Next time, I will introduce a program that applies displacement boundary conditions and calculates the overall stiffness equation (overall stiffness matrix).

Reference book

Shigeru Nagagi (2010) "First finite element method by Java" Corona Publishing Co., Ltd.

Recommended Posts

Finite element method Beginning 2 (Preparation of one-dimensional overall stiffness matrix)
Finite element method beginning (one-dimensional element stiffness matrix)
Memorandum of introduction of EXODUS, a data model of the finite element method (FEM)