FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
m_common_struct.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!-------------------------------------------------------------------------------
7 implicit none
8 integer, parameter, private :: kreal = kind(0.0d0)
9
11 character(len=128) :: sys_name
12 integer :: sys_type
15 integer :: node_id(3)
16 real(kind=kreal) :: coordsys(3, 3)
17 end type tlocalcoordsys
18
19 type(tlocalcoordsys), pointer, save :: g_localcoordsys(:) => null()
20
22 integer :: n_rot
23 type(trotcond), pointer :: conds(:)
24 end type
25
27 logical :: active
28 integer :: center_ngrp_id
29 integer :: torque_ngrp_id
30 real(kind=kreal) :: vec(3)
31 end type
32
33contains
34
37 integer :: i
38
39 if( .not. associated(g_localcoordsys) ) return
40 do i=1,size(g_localcoordsys)
41 g_localcoordsys(i)%sys_name = ""
42 g_localcoordsys(i)%sys_type = 10
43 g_localcoordsys(i)%node_ID(:) = -1
44 g_localcoordsys(i)%CoordSys(:,:) = 0.d0
45 enddo
46 end subroutine
47
49 subroutine print_localcoordsys(nfile, coordsys)
50 integer, intent(in) :: nfile
51 type(tlocalcoordsys), intent(in) :: coordsys
52
53 write(nfile, *) coordsys%sys_type, coordsys%sys_name
54 write(nfile, *) coordsys%node_ID(:)
55 write(nfile, *) coordsys%CoordSys(1,:)
56 write(nfile, *) coordsys%CoordSys(2,:)
57 write(nfile, *) coordsys%CoordSys(3,:)
58 end subroutine
59
61 logical function iscoordneeds( coordsys )
62 type(tlocalcoordsys), intent(in) :: coordsys
63
64 integer stype
65 iscoordneeds = .false.
66 stype = mod( coordsys%sys_type, 10 )
67 if( stype==0 ) iscoordneeds = .true.
68 end function
69
71 subroutine set_localcoordsys(coords, coordsys, outsys, ierr)
72 use m_utilities, only: cross_product
73 real(kind=kreal), intent(inout) :: coords(:, :) ! coord needs to define local coord sys
74 type(tlocalcoordsys), intent(in) :: coordsys
75 real(kind=kreal), intent(out) :: outsys(3, 3)
76 integer, intent(out) :: ierr
77
78 real(kind=kreal) :: f, xyza(3), xyzb(3), xyzc(3)
79 integer :: stype
80
81 ierr = 0
82 stype = mod( coordsys%sys_type, 10 )
83
84 if( stype==0 ) then ! given coordinates
85 outsys = coordsys%CoordSys
86 elseif( stype==1 ) then ! given nodes id
87 xyza = coords(1,:)-coords(3,:)
88 xyzb = coords(2,:)-coords(3,:)
89 call cross_product(xyza,xyzb,xyzc)
90 f = dsqrt( dot_product( xyza, xyza ) )
91 outsys(1,:) = xyza/f
92 f = dsqrt( dot_product( xyzc, xyzc ) )
93 outsys(3,:) = xyzc/f
94 call cross_product(outsys(3,:), outsys(1,:), outsys(2,:) )
95 endif
96 end subroutine
97
98 subroutine fstr_rotinfo_init(n, rinfo)
99 integer, intent(in) :: n
100 type(trotinfo), intent(inout) :: rinfo
101
102 integer :: i
103
104 if( n < 1 ) return
105
106 rinfo%n_rot = n
107 allocate(rinfo%conds(n))
108 do i=1,n
109 rinfo%conds(i)%active = .false.
110 rinfo%conds(i)%center_ngrp_id = -1
111 rinfo%conds(i)%torque_ngrp_id = -1
112 rinfo%conds(i)%vec(1:3) = 0.d0
113 enddo
114 end subroutine
115
116 subroutine fstr_rotinfo_finalize(rinfo)
117 type(trotinfo), intent(inout) :: rinfo
118
119 rinfo%n_rot = 0
120 if( associated(rinfo%conds)) deallocate(rinfo%conds)
121 end subroutine
122
123end module
This modules defines common structures for fem analysis.
logical function iscoordneeds(coordsys)
if need to fetch global nodes' coodinate
subroutine fstr_rotinfo_init(n, rinfo)
type(tlocalcoordsys), dimension(:), pointer, save g_localcoordsys
subroutine set_localcoordsys(coords, coordsys, outsys, ierr)
setup of coordinate system
subroutine fstr_rotinfo_finalize(rinfo)
subroutine print_localcoordsys(nfile, coordsys)
output of coordinate system
subroutine fstr_localcoordsys_init()
Initializer of global data.
This module provides aux functions.
Definition: utilities.f90:6
subroutine cross_product(v1, v2, vn)
Definition: utilities.f90:330