FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
hecmw_ML_helper.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!-------------------------------------------------------------------------------
5
6subroutine hecmw_ml_get_nlocal(id, nlocal, nlocal_allcolumns, ierr)
7 use hecmw_util
9 implicit none
10 integer(kind=kint), intent(in) :: id
11 integer(kind=kint), intent(out) :: nlocal
12 integer(kind=kint), intent(out) :: nlocal_allcolumns
13 integer(kind=kint), intent(out) :: ierr
14 type(hecmwst_matrix), pointer :: hecMAT
15 type(hecmwst_local_mesh), pointer :: hecMESH
16 call hecmw_mat_id_get(id, hecmat, hecmesh)
17 nlocal = hecmat%N * hecmat%NDOF
18 nlocal_allcolumns = hecmat%NP * hecmat%NDOF
19 ierr = 0
20end subroutine hecmw_ml_get_nlocal
21
22subroutine hecmw_ml_get_coord(id, x, y, z, ierr)
23 use hecmw_util
24 use hecmw_mat_id
25 implicit none
26 integer(kind=kint), intent(in) :: id
27 real(kind=kreal), intent(out) :: x(*), y(*), z(*)
28 integer(kind=kint), intent(out) :: ierr
29 type(hecmwst_matrix), pointer :: hecMAT
30 type(hecmwst_local_mesh), pointer :: hecMESH
31 integer(kind=kint) :: offset, i
32 call hecmw_mat_id_get(id, hecmat, hecmesh)
33 offset = 0
34 do i = 1, hecmesh%nn_internal
35 x(i) = hecmesh%node(offset+1)
36 y(i) = hecmesh%node(offset+2)
37 z(i) = hecmesh%node(offset+3)
38 offset = offset + 3
39 enddo
40 ierr = 0
41end subroutine hecmw_ml_get_coord
42
43subroutine hecmw_ml_get_rbm(id, rbm, ierr)
44 use hecmw_util
45 use hecmw_mat_id
46 use hecmw_etype
47 implicit none
48 integer(kind=kint), intent(in) :: id
49 real(kind=kreal), intent(out) :: rbm(*)
50 integer(kind=kint), intent(out) :: ierr
51 type(hecmwst_matrix), pointer :: hecMAT
52 type(hecmwst_local_mesh), pointer :: hecMESH
53 integer(kind=kint) :: Ndof, vec_leng, node, offset
54 real(kind=kreal) :: x, y, z
55
56 integer(kind=kint), allocatable :: mark(:)
57 integer(kind=kint) :: itype, ic_type, nn, is, iE, icel, iiS, j, nod
58
59 call hecmw_mat_id_get(id, hecmat, hecmesh)
60
61 ! Mark nodes for rotational DOF
62 allocate(mark(hecmesh%n_node))
63 mark = 0
64 do itype = 1, hecmesh%n_elem_type
65 ic_type = hecmesh%elem_type_item(itype)
66 if (hecmw_is_etype_33struct(ic_type)) then
67 nn = hecmw_get_max_node(ic_type)
68 is = hecmesh%elem_type_index(itype-1)+1
69 ie = hecmesh%elem_type_index(itype )
70 do icel = is, ie
71 iis = hecmesh%elem_node_index(icel-1)
72 ! mark latter halves of the nodes
73 do j = nn/2+1, nn
74 nod = hecmesh%elem_node_item(iis+j)
75 mark(nod) = 1
76 enddo
77 enddo
78 endif
79 enddo
80 ndof = hecmat%NDOF
81 vec_leng = hecmesh%nn_internal * ndof
82
83 if (ndof == 1) then
84
85 rbm(1:vec_leng)=1.d0
86
87 else if (ndof == 2) then
88
89 rbm(1:vec_leng*3)=0.0d0
90
91 do node = 1, hecmesh%nn_internal
92 x = hecmesh%node(3*node-2)
93 y = hecmesh%node(3*node-1)
94
95 ! translation x
96 offset = (node-1)*ndof
97 rbm(offset+1)=1.d0
98 rbm(offset+2)=0.d0
99
100 ! translation y
101 offset = offset + vec_leng
102 rbm(offset+1)=0.d0
103 rbm(offset+2)=1.d0
104
105 ! rotation z
106 offset = offset + vec_leng
107 rbm(offset+1)= -y
108 rbm(offset+2)= x
109
110 enddo
111
112 else
113
114 rbm(1:vec_leng*6)=0.0d0
115 do node = 1, hecmesh%nn_internal
116 if (mark(node) == 0) then
117!!! translational DOF
118
119 x = hecmesh%node(3*node-2)
120 y = hecmesh%node(3*node-1)
121 z = hecmesh%node(3*node )
122
123 ! translation x
124 offset = (node-1)*ndof
125 rbm(offset+1)=1.d0
126 rbm(offset+2)=0.d0
127 rbm(offset+3)=0.d0
128
129 ! translation y
130 offset = offset + vec_leng
131 rbm(offset+1)=0.d0
132 rbm(offset+2)=1.d0
133 rbm(offset+3)=0.d0
134
135 ! translation z
136 offset = offset + vec_leng
137 rbm(offset+1)=0.d0
138 rbm(offset+2)=0.d0
139 rbm(offset+3)=1.d0
140
141 ! rotation x
142 offset = offset + vec_leng
143 rbm(offset+1)=0.d0
144 rbm(offset+2)= -z
145 rbm(offset+3)= y
146
147 ! rotation y
148 offset = offset + vec_leng
149 rbm(offset+1)= z
150 rbm(offset+2)=0.d0
151 rbm(offset+3)= -x
152
153 ! rotation z
154 offset = offset + vec_leng
155 rbm(offset+1)= -y
156 rbm(offset+2)= x
157 rbm(offset+3)=0.d0
158
159 else
160!!! rotational DOF
161
162 ! translation x
163 offset = (node-1)*ndof
164 rbm(offset+1)=0.d0
165 rbm(offset+2)=0.d0
166 rbm(offset+3)=0.d0
167
168 ! translation y
169 offset = offset + vec_leng
170 rbm(offset+1)=0.d0
171 rbm(offset+2)=0.d0
172 rbm(offset+3)=0.d0
173
174 ! translation z
175 offset = offset + vec_leng
176 rbm(offset+1)=0.d0
177 rbm(offset+2)=0.d0
178 rbm(offset+3)=0.d0
179
180 ! rotation x
181 offset = offset + vec_leng
182 rbm(offset+1)=1.d0
183 rbm(offset+2)=0.d0
184 rbm(offset+3)=0.d0
185
186 ! rotation y
187 offset = offset + vec_leng
188 rbm(offset+1)=0.d0
189 rbm(offset+2)=1.d0
190 rbm(offset+3)=0.d0
191
192 ! rotation z
193 offset = offset + vec_leng
194 rbm(offset+1)=0.d0
195 rbm(offset+2)=0.d0
196 rbm(offset+3)=1.d0
197 endif
198 enddo
199 endif
200
201 deallocate(mark)
202 ierr = 0
203end subroutine hecmw_ml_get_rbm
204
205subroutine hecmw_ml_get_loglevel(id, level)
206 use hecmw_util
208 use hecmw_mat_id
209 implicit none
210 integer(kind=kint), intent(in) :: id
211 integer(kind=kint), intent(out) :: level
212 type(hecmwst_matrix), pointer :: hecMAT
213 type(hecmwst_local_mesh), pointer :: hecMESH
214 call hecmw_mat_id_get(id, hecmat, hecmesh)
215 level = hecmw_mat_get_timelog(hecmat)
216end subroutine hecmw_ml_get_loglevel
217
218subroutine hecmw_ml_get_opt1(id, opt1, ierr)
219 use hecmw_util
220 use hecmw_mat_id
222 implicit none
223 integer(kind=kint), intent(in) :: id
224 integer(kind=kint), intent(out) :: opt1
225 integer(kind=kint), intent(out) :: ierr
226 type(hecmwst_matrix), pointer :: hecMAT
227 type(hecmwst_local_mesh), pointer :: hecMESH
228 call hecmw_mat_id_get(id, hecmat, hecmesh)
229 opt1 = hecmw_mat_get_solver_opt1(hecmat)
230 ierr = 0
231end subroutine hecmw_ml_get_opt1
232
233subroutine hecmw_ml_get_opt2(id, opt2, ierr)
234 use hecmw_util
235 use hecmw_mat_id
237 implicit none
238 integer(kind=kint), intent(in) :: id
239 integer(kind=kint), intent(out) :: opt2
240 integer(kind=kint), intent(out) :: ierr
241 type(hecmwst_matrix), pointer :: hecMAT
242 type(hecmwst_local_mesh), pointer :: hecMESH
243 call hecmw_mat_id_get(id, hecmat, hecmesh)
244 opt2 = hecmw_mat_get_solver_opt2(hecmat)
245 ierr = 0
246end subroutine hecmw_ml_get_opt2
247
248subroutine hecmw_ml_get_opt3(id, opt3, ierr)
249 use hecmw_util
250 use hecmw_mat_id
252 implicit none
253 integer(kind=kint), intent(in) :: id
254 integer(kind=kint), intent(out) :: opt3
255 integer(kind=kint), intent(out) :: ierr
256 type(hecmwst_matrix), pointer :: hecMAT
257 type(hecmwst_local_mesh), pointer :: hecMESH
258 call hecmw_mat_id_get(id, hecmat, hecmesh)
259 opt3 = hecmw_mat_get_solver_opt3(hecmat)
260 ierr = 0
261end subroutine hecmw_ml_get_opt3
262
263subroutine hecmw_ml_get_opt4(id, opt4, ierr)
264 use hecmw_util
265 use hecmw_mat_id
267 implicit none
268 integer(kind=kint), intent(in) :: id
269 integer(kind=kint), intent(out) :: opt4
270 integer(kind=kint), intent(out) :: ierr
271 type(hecmwst_matrix), pointer :: hecMAT
272 type(hecmwst_local_mesh), pointer :: hecMESH
273 call hecmw_mat_id_get(id, hecmat, hecmesh)
274 opt4 = hecmw_mat_get_solver_opt4(hecmat)
275 ierr = 0
276end subroutine hecmw_ml_get_opt4
277
278subroutine hecmw_ml_get_opt5(id, opt5, ierr)
279 use hecmw_util
280 use hecmw_mat_id
282 implicit none
283 integer(kind=kint), intent(in) :: id
284 integer(kind=kint), intent(out) :: opt5
285 integer(kind=kint), intent(out) :: ierr
286 type(hecmwst_matrix), pointer :: hecMAT
287 type(hecmwst_local_mesh), pointer :: hecMESH
288 call hecmw_mat_id_get(id, hecmat, hecmesh)
289 opt5 = hecmw_mat_get_solver_opt5(hecmat)
290 ierr = 0
291end subroutine hecmw_ml_get_opt5
292
293subroutine hecmw_ml_get_opt6(id, opt6, ierr)
294 use hecmw_util
295 use hecmw_mat_id
297 implicit none
298 integer(kind=kint), intent(in) :: id
299 integer(kind=kint), intent(out) :: opt6
300 integer(kind=kint), intent(out) :: ierr
301 type(hecmwst_matrix), pointer :: hecMAT
302 type(hecmwst_local_mesh), pointer :: hecMESH
303 call hecmw_mat_id_get(id, hecmat, hecmesh)
304 opt6 = hecmw_mat_get_solver_opt6(hecmat)
305 ierr = 0
306end subroutine hecmw_ml_get_opt6
307
308subroutine hecmw_ml_set_opt6(id, opt6, ierr)
309 use hecmw_util
310 use hecmw_mat_id
312 implicit none
313 integer(kind=kint), intent(in) :: id
314 integer(kind=kint), intent(in) :: opt6
315 integer(kind=kint), intent(out) :: ierr
316 type(hecmwst_matrix), pointer :: hecMAT
317 type(hecmwst_local_mesh), pointer :: hecMESH
318 call hecmw_mat_id_get(id, hecmat, hecmesh)
319 call hecmw_mat_set_solver_opt6(hecmat, opt6)
320 ierr = 0
321end subroutine hecmw_ml_set_opt6
subroutine hecmw_ml_set_opt6(id, opt6, ierr)
subroutine hecmw_ml_get_rbm(id, rbm, ierr)
subroutine hecmw_ml_get_opt3(id, opt3, ierr)
subroutine hecmw_ml_get_opt1(id, opt1, ierr)
subroutine hecmw_ml_get_opt2(id, opt2, ierr)
subroutine hecmw_ml_get_opt6(id, opt6, ierr)
subroutine hecmw_ml_get_opt5(id, opt5, ierr)
subroutine hecmw_ml_get_coord(id, x, y, z, ierr)
subroutine hecmw_ml_get_nlocal(id, nlocal, nlocal_allcolumns, ierr)
subroutine hecmw_ml_get_opt4(id, opt4, ierr)
subroutine hecmw_ml_get_loglevel(id, level)
I/O and Utility.
integer(kind=kint) function hecmw_get_max_node(etype)
logical function hecmw_is_etype_33struct(etype)
subroutine, public hecmw_mat_id_get(id, hecmat, hecmesh)
integer(kind=kint) function, public hecmw_mat_get_solver_opt6(hecmat)
integer(kind=kint) function, public hecmw_mat_get_solver_opt4(hecmat)
integer(kind=kint) function, public hecmw_mat_get_solver_opt2(hecmat)
subroutine, public hecmw_mat_set_solver_opt6(hecmat, solver_opt6)
integer(kind=kint) function, public hecmw_mat_get_timelog(hecmat)
integer(kind=kint) function, public hecmw_mat_get_solver_opt3(hecmat)
integer(kind=kint) function, public hecmw_mat_get_solver_opt1(hecmat)
integer(kind=kint) function, public hecmw_mat_get_solver_opt5(hecmat)
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal