Coverage for src / sdynpy / demo / beam_plate.py: 100%

1# -*- coding: utf-8 -*- 

2""" 

3Demonstration system for a plate-like structure made of beams 

4""" 

5""" 

6Copyright 2022 National Technology & Engineering Solutions of Sandia, 

7LLC (NTESS). Under the terms of Contract DE-NA0003525 with NTESS, the U.S. 

8Government retains certain rights in this software. 

9 

10This program is free software: you can redistribute it and/or modify 

11it under the terms of the GNU General Public License as published by 

12the Free Software Foundation, either version 3 of the License, or 

13(at your option) any later version. 

14 

15This program is distributed in the hope that it will be useful, 

16but WITHOUT ANY WARRANTY; without even the implied warranty of 

17MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 

18GNU General Public License for more details. 

19 

20You should have received a copy of the GNU General Public License 

21along with this program. If not, see <https://www.gnu.org/licenses/>. 

22""" 

23 

24import numpy as _np 

25from ..fem.sdynpy_beam import beamkm as _beamkm 

26from ..core.sdynpy_geometry import (Geometry as _Geometry, 

27 node_array as _node_array, 

28 traceline_array as _traceline_array, 

29 coordinate_system_array as _coordinate_system_array) 

30from ..core.sdynpy_coordinate import from_nodelist as _from_nodelist 

31from ..core.sdynpy_system import (System as _System, 

32 substructure_by_position as _substructure_by_position) 

33 

34 

35def create_models(width=1, length=0.5, grid=0.125, 

36 beam_height=0.02, beam_width=0.02, 

37 E=69e9, nu=0.33, rho=2830): 

38 nodes_width = int(_np.round(width / grid)) + 1 

39 nodes_length = int(_np.round(length / grid)) + 1 

40 coords = _np.moveaxis(_np.meshgrid(_np.linspace(0, width, nodes_width), 

41 _np.linspace(0, length, nodes_length), 

42 0, 

43 indexing='ij'), 0, -1)[:, :, 0, :] 

44 indices = _np.arange(_np.prod(coords.shape[:-1])).reshape(coords.shape[:-1]) 

45 

46 conn = [] 

47 for i in range(coords.shape[0]): 

48 for j in range(coords.shape[1]): 

49 try: 

50 conn.append([indices[i, j], indices[i, j + 1]]) 

51 except IndexError: 

52 pass 

53 try: 

54 conn.append([indices[i, j], indices[i + 1, j]]) 

55 except IndexError: 

56 pass 

57 conn = _np.array(conn) 

58 

59 bend_direction_1 = _np.array( 

60 [_np.zeros(conn.shape[0]), 

61 _np.zeros(conn.shape[0]), 

62 _np.ones(conn.shape[0])] 

63 ).T 

64 

65 I1 = beam_width * beam_height**3 / 12 

66 I2 = beam_width**3 * beam_height / 12 

67 G = E / (2 * (1 - nu)) 

68 J = I1 + I2 

69 Ixx_per_L = (1 / 12) * rho * beam_width * beam_height * (beam_width**2 + beam_height**2) 

70 

71 # Create arguments for the beamkm function 

72 ae = beam_width * beam_height * E * _np.ones(conn.shape[0]) 

73 jg = J * G * _np.ones(conn.shape[0]) 

74 ei1 = E * I1 * _np.ones(conn.shape[0]) 

75 ei2 = E * I2 * _np.ones(conn.shape[0]) 

76 mass_per_length = rho * beam_width * beam_height * _np.ones(conn.shape[0]) 

77 tmmi_per_length = Ixx_per_L * _np.ones(conn.shape[0]) 

78 

79 # Compute K and M via beamkm 

80 this_K, this_M = _beamkm(coords.reshape(-1, 3), conn, bend_direction_1, 

81 ae, jg, ei1, ei2, mass_per_length, tmmi_per_length) 

82 

83 geometry = _Geometry( 

84 _node_array(_np.arange(coords.reshape(-1, 3).shape[0]) + 1, coords.reshape(-1, 3), 1), 

85 _coordinate_system_array(), 

86 _traceline_array(_np.arange(conn.shape[0]) + 1, color=1, connectivity=conn + 1) 

87 ) 

88 

89 system = _System( 

90 _from_nodelist(geometry.node.id, [1, 2, 3, 4, 5, 6]), 

91 this_M, this_K 

92 ) 

93 

94 return system, geometry 

95 

96 

97system, geometry = create_models()