FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
static_echo.f90
Go to the documentation of this file.
1!-------------------------------------------------------------------------------
2! Copyright (c) 2019 FrontISTR Commons
3! This software is released under the MIT License, see LICENSE.txt
4!-------------------------------------------------------------------------------
7contains
8 !C
9 !C***
11 !C***
12 !C
13 subroutine fstr_echo (hecMESH)
14 use m_fstr
16 implicit none
17 type (hecmwST_local_mesh):: hecMESH
18
19 !** Local variables
20 integer(kind=kint) :: i,j,is,iE,nn,ic_type
21 integer(kind=kint) :: icel,isect,ig1,iS0,iE0,ik
22 integer(kind=kint) :: in,nid, itype
23 real(kind=kreal) :: x,y,z
24 integer(kind=kint):: nids(20)
25 !C====
26 !C +-------------------------------+
27 !C | NODE |
28 !C +-------------------------------+
29 !C===
30 !C** nodal coordinate
31 write(ilog,*) '### Number of nodes',hecmesh%n_node
32 write(ilog,*) 'ID X Y Z'
33 do i=1,hecmesh%n_node
34 nid = hecmesh%global_node_ID(i)
35 x = hecmesh%node(3*i-2)
36 y = hecmesh%node(3*i-1)
37 z = hecmesh%node(3*i)
38 write(ilog,'(i8,e15.5,e15.5,e15.5)') nid,x,y,z
39 enddo
40 !C
41 !C +-------------------------------+
42 !C | ELEMENT |
43 !C +-------------------------------+
44 !C===
45 call fstr2hecmw_mesh_conv( hecmesh )
46 write(ilog,*) '### Elements', hecmesh%n_elem
47 do itype= 1, hecmesh%n_elem_type
48 is= hecmesh%elem_type_index(itype-1) + 1
49 ie= hecmesh%elem_type_index(itype )
50 ic_type= hecmesh%elem_type_item(itype)
51 !C** Set number of nodes
52 nn = hecmw_get_max_node(ic_type)
53 !C** element loop
54 do icel= is, ie
55 !C** node ID
56 is= hecmesh%elem_node_index(icel-1)
57 do j=1,nn
58 if( hecmesh%n_refine > 0 ) then
59 nids(j)= hecmesh%elem_node_item (is+j)
60 else
61 nids(j)= hecmesh%global_node_ID( hecmesh%elem_node_item (is+j))
62 endif
63 enddo
64 !C** section ID
65 isect= hecmesh%section_ID(icel)
66 write(ilog,*) '### Element ID=',ic_type,isect,hecmesh%global_elem_id(icel)
67 write(ilog,*) (nids(j),j=1,nn)
68 enddo
69 enddo
70 call hecmw2fstr_mesh_conv( hecmesh )
71 !C +-------------------------------+
72 !C | NODE GROUP |
73 !C +-------------------------------+
74 write(ilog,*) '### Ngroup'
75 do ig1= 1, hecmesh%node_group%n_grp
76 write(ilog,*)
77 write(ilog,'(a80)') hecmesh%node_group%grp_name(ig1)
78 is0= hecmesh%node_group%grp_index(ig1-1) + 1
79 ie0= hecmesh%node_group%grp_index(ig1 )
80 do ik= is0, ie0
81 in = hecmesh%node_group%grp_item(ik)
82 write(ilog,*) hecmesh%global_node_ID(in)
83 enddo
84 enddo
85 !C +-------------------------------+
86 !C | ELEMEN GROUP |
87 !C +-------------------------------+
88 write(ilog,*) '### Egroup'
89 do ig1= 1, hecmesh%elem_group%n_grp
90 write(ilog,*)
91 write(ilog,'(a80)') hecmesh%elem_group%grp_name(ig1)
92 is0= hecmesh%elem_group%grp_index(ig1-1) + 1
93 ie0= hecmesh%elem_group%grp_index(ig1 )
94 do ik= is0, ie0
95 in = hecmesh%elem_group%grp_item(ik)
96 write(ilog,*) hecmesh%global_elem_ID(in)
97 enddo
98 enddo
99 write(ilog,*) '### Reftemp',ref_temp
100 !C====
101 end subroutine fstr_echo
102end module m_static_echo
This module defined coomon data and basic structures for analysis.
Definition: m_fstr.f90:15
integer(kind=kint), parameter ilog
FILE HANDLER.
Definition: m_fstr.f90:91
real(kind=kreal), pointer ref_temp
REFTEMP.
Definition: m_fstr.f90:120
HECMW to FSTR Mesh Data Converter. Convering Conectivity of Element Type 232, 342 and 352.
subroutine hecmw2fstr_mesh_conv(hecmesh)
subroutine fstr2hecmw_mesh_conv(hecmesh)
This module provide a function to ECHO for IFSTR solver.
Definition: static_echo.f90:6
subroutine fstr_echo(hecmesh)
ECHO for IFSTR solver.
Definition: static_echo.f90:14