FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
fstr_mat_con_contact.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!-------------------------------------------------------------------------------
6! stiffness matrix structure for the contact analysis
7! employing standard Lagrange multiplier algorithm
8
10
11 use m_fstr
12 use elementinfo
13 use mcontact
14
15 implicit none
16
19 integer(kind=kint) :: num_node = 0, num_lagrange = 0
20 integer(kind=kint), pointer :: id_node(:) => null()
21 integer(kind=kint), pointer :: id_lagrange(:) => null()
22 end type
23
27 integer(kind=kint) :: num_lagrange = 0
28 integer(kind=kint) :: numl_lagrange = 0, numu_lagrange = 0
30 integer(kind=kint), pointer :: indexl_lagrange(:) => null(), &
31 indexu_lagrange(:) => null()
32 integer(kind=kint), pointer :: iteml_lagrange(:) => null(), &
33 itemu_lagrange(:) => null()
34 real(kind=kreal), pointer :: al_lagrange(:) => null(), &
35 au_lagrange(:) => null()
36 real(kind=kreal), pointer :: lagrange(:) => null()
38
39 integer(kind=kint), save :: npl_org, npu_org
40 type(noderelated), pointer, save :: list_noderelated_org(:) => null()
41
42 type(noderelated), pointer :: list_noderelated(:) => null()
43
44 logical :: permission = .false.
45
46 private :: insert_lagrange, insert_node, find_locationinarray, reallocate_memory
47
48contains
49
50
53
54 type(hecmwst_matrix) :: hecMAT
55
56 integer(kind=kint) :: numL, numU, num_nodeRelated
57 integer(kind=kint) :: i, j
58 integer(kind=kint) :: ierr
59
60 npl_org = hecmat%NPL; npu_org = hecmat%NPU
61
62 if( associated(list_noderelated_org) ) return
63 allocate(list_noderelated_org(hecmat%NP), stat=ierr)
64 if( ierr /= 0) stop " Allocation error, list_nodeRelated_org "
65
66 do i = 1, hecmat%NP
67
68 numl = hecmat%indexL(i) - hecmat%indexL(i-1)
69 numu = hecmat%indexU(i) - hecmat%indexU(i-1)
70
71 num_noderelated = numl + numu + 1
72
73 allocate(list_noderelated_org(i)%id_node(num_noderelated), stat=ierr)
74 if( ierr /= 0) stop " Allocation error, list_nodeRelated_org%id_node "
75
76 list_noderelated_org(i)%num_node = num_noderelated
77
78 do j = 1, numl
79 list_noderelated_org(i)%id_node(j) = hecmat%itemL(hecmat%indexL(i-1)+j)
80 enddo
81 list_noderelated_org(i)%id_node(numl+1) = i
82 do j = 1, numu
83 list_noderelated_org(i)%id_node(numl+1+j) = hecmat%itemU(hecmat%indexU(i-1)+j)
84 enddo
85
86 enddo
87
89
92 subroutine fstr_mat_con_contact(cstep,hecMAT,fstrSOLID,fstrMAT,infoCTChange,conMAT)
93
94 integer(kind=kint) :: cstep
95 type(hecmwst_matrix) :: hecMAT
96 type(fstr_solid) :: fstrSOLID
97 type(fstrst_matrix_contact_lagrange) :: fstrMAT
98 type(fstr_info_contactchange) :: infoCTChange
99
100 integer(kind=kint) :: num_lagrange
101 integer(kind=kint) :: countNon0LU_node, countNon0LU_lagrange
102 integer(kind=kint) :: numNon0_node, numNon0_lagrange
104 type (hecmwST_matrix),optional :: conMAT
105
106 num_lagrange = infoctchange%contactNode_current
107 fstrmat%num_lagrange = num_lagrange
108
109 ! Get original list of related nodes
110 call getoriginallistofrelatednodes(hecmat%NP,num_lagrange)
111
112 ! Construct new list of related nodes and Lagrange multipliers
113 countnon0lu_node = npl_org + npu_org
114 countnon0lu_lagrange = 0
115 if( fstr_is_contact_active() ) &
116 call getnewlistofrelatednodesandlagrangemultipliers(cstep,hecmat%NP,fstrsolid,countnon0lu_node,countnon0lu_lagrange)
117
118 ! Construct new matrix structure(hecMAT&fstrMAT)
119 numnon0_node = countnon0lu_node/2
120 numnon0_lagrange = countnon0lu_lagrange/2
121 ! ---- For Parallel Contact with Multi-Partition Domains
122 if(paracontactflag.and.present(conmat)) then
123 call constructnewmatrixstructure(hecmat,fstrmat,numnon0_node,numnon0_lagrange,conmat)
124 else
125 call constructnewmatrixstructure(hecmat,fstrmat,numnon0_node,numnon0_lagrange)
126 endif
127
128 ! Copy Lagrange multipliers
129 if( fstr_is_contact_active() ) &
130 call fstr_copy_lagrange_contact(fstrsolid,fstrmat)
131
132 end subroutine fstr_mat_con_contact
133
135 subroutine getoriginallistofrelatednodes(np,num_lagrange)
136
137 integer(kind=kint) :: np, num_lagrange
138 integer(kind=kint) :: num_nodeRelated_org
139 integer(kind=kint) :: i, ierr
140
141 allocate(list_noderelated(np+num_lagrange), stat=ierr)
142 if( ierr /= 0) stop " Allocation error, list_nodeRelated "
143
144 do i = 1, np !hecMAT%NP
145 num_noderelated_org = list_noderelated_org(i)%num_node
146 allocate(list_noderelated(i)%id_node(num_noderelated_org), stat=ierr)
147 if( ierr /= 0) stop " Allocation error, list_nodeRelated%id_node "
148 list_noderelated(i)%num_node = num_noderelated_org
149 list_noderelated(i)%id_node(1:num_noderelated_org) = list_noderelated_org(i)%id_node(1:num_noderelated_org)
150 enddo
151
152 if( fstr_is_contact_active() ) then
153 do i = np+1, np+num_lagrange !hecMAT%NP+1, hecMAT%NP+num_lagrange
154 allocate(list_noderelated(i)%id_lagrange(5), stat=ierr)
155 if( ierr /= 0) stop " Allocation error, list_nodeRelated%id_lagrange "
156 list_noderelated(i)%num_lagrange = 0
157 list_noderelated(i)%id_lagrange = 0
158 enddo
159 endif
160
161 end subroutine getoriginallistofrelatednodes
162
163
165 subroutine getnewlistofrelatednodesandlagrangemultipliers(cstep,np,fstrSOLID,countNon0LU_node,countNon0LU_lagrange)
166
167 type(fstr_solid) :: fstrSOLID
168
169 integer(kind=kint) :: cstep
170 integer(kind=kint) :: np
171 integer(kind=kint) :: countNon0LU_node, countNon0LU_lagrange
172 integer(kind=kint) :: grpid
173 integer(kind=kint) :: count_lagrange
174 integer(kind=kint) :: ctsurf, etype, nnode, ndLocal(l_max_surface_node + 1)
175 integer(kind=kint) :: i, j, k, l, num, num_nodeRelated_org, ierr
176 real(kind=kreal) :: fcoeff
177
178 count_lagrange = 0
179 do i = 1, size(fstrsolid%contacts)
180
181 grpid = fstrsolid%contacts(i)%group
182 if( .not. fstr_iscontactactive( fstrsolid, grpid, cstep ) ) cycle
183
184 fcoeff = fstrsolid%contacts(i)%fcoeff
185
186 do j = 1, size(fstrsolid%contacts(i)%slave)
187
188 if( fstrsolid%contacts(i)%states(j)%state == contactfree ) cycle
189 ctsurf = fstrsolid%contacts(i)%states(j)%surface
190 etype = fstrsolid%contacts(i)%master(ctsurf)%etype
191 if( etype/=fe_tri3n .and. etype/=fe_quad4n ) &
192 stop " ##Error: This element type is not supported in contact analysis !!! "
193 nnode = size(fstrsolid%contacts(i)%master(ctsurf)%nodes)
194 ndlocal(1) = fstrsolid%contacts(i)%slave(j)
195 ndlocal(2:nnode+1) = fstrsolid%contacts(i)%master(ctsurf)%nodes(1:nnode)
196
197 count_lagrange = count_lagrange + 1
198
199 do k = 1, nnode+1
200
201 if( .not. associated(list_noderelated(ndlocal(k))%id_lagrange) )then
202 num = 10
203 ! if( k == 1 ) num = 1
204 allocate(list_noderelated(ndlocal(k))%id_lagrange(num),stat=ierr)
205 if( ierr /= 0) stop " Allocation error, list_nodeRelated%id_lagrange "
206 list_noderelated(ndlocal(k))%num_lagrange = 0
207 list_noderelated(ndlocal(k))%id_lagrange = 0
208 endif
209
210 if( fcoeff /= 0.0d0 ) then
211 num_noderelated_org = list_noderelated_org(ndlocal(k))%num_node
212 if( list_noderelated(ndlocal(k))%num_node == num_noderelated_org )then
213 num = 10
214 if(k==1) num = 4
215 call reallocate_memory(num,list_noderelated(ndlocal(k)))
216 endif
217 endif
218
219 call insert_lagrange(k,count_lagrange,list_noderelated(ndlocal(k)),countnon0lu_lagrange)
220
221 do l = k, nnode+1
222 if( fcoeff /= 0.0d0 ) then
223 if( k /= l) then
224 num_noderelated_org = list_noderelated_org(ndlocal(k))%num_node
225 call insert_node(ndlocal(l),list_noderelated(ndlocal(k)),countnon0lu_node)
226 num_noderelated_org = list_noderelated_org(ndlocal(l))%num_node
227 if( list_noderelated(ndlocal(l))%num_node == num_noderelated_org )then
228 num = 10
229 call reallocate_memory(num,list_noderelated(ndlocal(l)))
230 endif
231 call insert_node(ndlocal(k),list_noderelated(ndlocal(l)),countnon0lu_node)
232 endif
233 endif
234
235 if(k == 1) &
236 call insert_lagrange(0,ndlocal(l),list_noderelated(np+count_lagrange),countnon0lu_lagrange)
237
238 enddo
239
240 enddo
241
242 enddo
243
244 enddo
245
247
248
250 subroutine constructnewmatrixstructure(hecMAT,fstrMAT,numNon0_node,numNon0_lagrange,conMAT)
251
252 type(hecmwst_matrix) :: hecmat
253 type(fstrst_matrix_contact_lagrange) :: fstrMAT
254
255 integer(kind=kint) :: numNon0_node, numNon0_lagrange
256 integer(kind=kint) :: countNon0L_node, countNon0U_node, countNon0U_lagrange, countNon0L_lagrange
257 integer(kind=kint) :: i, j, ierr
258 integer(kind=kint) :: numI_node, numI_lagrange
259 integer(kind=kint) :: ndof, nn
260 type(hecmwst_matrix), optional :: conMAT
261
262 ! ---- For Parallel Contact with Multi-Partition Domains
263 if(paracontactflag.and.present(conmat)) then
264 conmat%N = hecmat%N
265 conmat%NP = hecmat%NP
266 conmat%ndof = hecmat%ndof
267 if(associated(conmat%indexL).and.associated(conmat%indexU))deallocate(conmat%indexL,conmat%indexU)
268 allocate(conmat%indexL(0:conmat%NP), conmat%indexU(0:conmat%NP), stat=ierr)
269 if ( ierr /= 0) stop " Allocation error, conMAT%indexL-conMAT%indexU "
270 conmat%indexL = 0 ; conmat%indexU = 0
271 if(associated(conmat%itemL).and.associated(conmat%itemU))deallocate(conmat%itemL,conmat%itemU)
272 allocate(conmat%itemL(numnon0_node), conmat%itemU(numnon0_node), stat=ierr)
273 if ( ierr /= 0) stop " Allocation error, conMAT%itemL-conMAT%itemU "
274 conmat%itemL = 0 ; conmat%itemU = 0
275 !
276 conmat%NPL = numnon0_node
277 conmat%NPU = numnon0_node
278 endif
279
280 if(associated(hecmat%indexL).and.associated(hecmat%indexU))deallocate(hecmat%indexL,hecmat%indexU)
281 allocate(hecmat%indexL(0:hecmat%NP), hecmat%indexU(0:hecmat%NP), stat=ierr)
282 if ( ierr /= 0) stop " Allocation error, hecMAT%indexL-hecMAT%indexU "
283 hecmat%indexL = 0 ; hecmat%indexU = 0
284 if(associated(hecmat%itemL).and.associated(hecmat%itemU))deallocate(hecmat%itemL,hecmat%itemU)
285 allocate(hecmat%itemL(numnon0_node), hecmat%itemU(numnon0_node), stat=ierr)
286 if ( ierr /= 0) stop " Allocation error, hecMAT%itemL-hecMAT%itemU "
287 hecmat%itemL = 0 ; hecmat%itemU = 0
288
289 if(associated(fstrmat%indexL_lagrange).and.associated(fstrmat%indexU_lagrange)) &
290 deallocate(fstrmat%indexL_lagrange,fstrmat%indexU_lagrange)
291 if(associated(fstrmat%itemL_lagrange).and.associated(fstrmat%itemU_lagrange)) &
292 deallocate(fstrmat%itemL_lagrange,fstrmat%itemU_lagrange)
293 if( fstr_is_contact_active() ) then
294 allocate(fstrmat%indexL_lagrange(0:fstrmat%num_lagrange), fstrmat%indexU_lagrange(0:hecmat%NP), stat=ierr)
295 if ( ierr /= 0) stop " Allocation error, fstrMAT%indexL_lagrange-fstrMAT%indexU_lagrange "
296 fstrmat%indexL_lagrange = 0 ; fstrmat%indexU_lagrange = 0
297 allocate(fstrmat%itemL_lagrange(numnon0_lagrange), fstrmat%itemU_lagrange(numnon0_lagrange), stat=ierr)
298 if ( ierr /= 0) stop " Allocation error, fstrMAT%itemL_lagrange-fstrMAT%itemU_lagrange "
299 fstrmat%itemL_lagrange = 0 ; fstrmat%itemU_lagrange = 0
300 endif
301
302 hecmat%NPL = numnon0_node
303 hecmat%NPU = numnon0_node
304
305 fstrmat%numL_lagrange = numnon0_lagrange
306 fstrmat%numU_lagrange = numnon0_lagrange
307
308 countnon0l_node = 0
309 countnon0u_node = 0
310 countnon0u_lagrange = 0
311 do i = 1, hecmat%NP
312
313 list_noderelated(i)%num_node = count(list_noderelated(i)%id_node /= 0)
314 numi_node = list_noderelated(i)%num_node
315 if( fstr_is_contact_active() ) &
316 numi_lagrange = list_noderelated(i)%num_lagrange
317
318 do j = 1, numi_node
319 if( list_noderelated(i)%id_node(j) < i ) then
320 countnon0l_node = countnon0l_node + 1
321 hecmat%itemL(countnon0l_node) = list_noderelated(i)%id_node(j)
322 elseif( list_noderelated(i)%id_node(j) > i ) then
323 countnon0u_node = countnon0u_node + 1
324 hecmat%itemU(countnon0u_node) = list_noderelated(i)%id_node(j)
325 endif
326 enddo
327 hecmat%indexL(i) = countnon0l_node
328 hecmat%indexU(i) = countnon0u_node
329
330 if( fstr_is_contact_active() ) then
331 do j = 1, numi_lagrange
332 countnon0u_lagrange = countnon0u_lagrange + 1
333 fstrmat%itemU_lagrange(countnon0u_lagrange) = list_noderelated(i)%id_lagrange(j)
334 enddo
335 fstrmat%indexU_lagrange(i) = countnon0u_lagrange
336 endif
337
338 deallocate(list_noderelated(i)%id_node)
339 if(associated(list_noderelated(i)%id_lagrange)) deallocate(list_noderelated(i)%id_lagrange)
340
341 end do
342
343 ! ---- For Parallel Contact with Multi-Partition Domains
344 if(paracontactflag.and.present(conmat)) then
345 conmat%itemL(:) = hecmat%itemL(:)
346 conmat%indexL(:) = hecmat%indexL(:)
347 conmat%itemU(:) = hecmat%itemU(:)
348 conmat%indexU(:) = hecmat%indexU(:)
349 endif
350
351 if( fstr_is_contact_active() ) then
352 countnon0l_lagrange = 0
353 do i = 1, fstrmat%num_lagrange
354 numi_lagrange = list_noderelated(hecmat%NP+i)%num_lagrange
355 do j = 1, numi_lagrange
356 countnon0l_lagrange = countnon0l_lagrange + 1
357 fstrmat%itemL_lagrange(countnon0l_lagrange) = list_noderelated(hecmat%NP+i)%id_lagrange(j)
358 enddo
359 fstrmat%indexL_lagrange(i) = countnon0l_lagrange
360 deallocate(list_noderelated(hecmat%NP+i)%id_lagrange)
361 enddo
362 endif
363
364 deallocate(list_noderelated)
365
366 ndof = hecmat%NDOF
367 nn = ndof*ndof
368 if(associated(hecmat%AL)) deallocate(hecmat%AL)
369 allocate(hecmat%AL(nn*hecmat%NPL), stat=ierr)
370 if ( ierr /= 0 ) stop " Allocation error, hecMAT%AL "
371 hecmat%AL = 0.0d0
372
373 if(associated(hecmat%AU)) deallocate(hecmat%AU)
374 allocate(hecmat%AU(nn*hecmat%NPU), stat=ierr)
375 if ( ierr /= 0 ) stop " Allocation error, hecMAT%AU "
376 hecmat%AU = 0.0d0
377
378 if(associated(fstrmat%AL_lagrange)) deallocate(fstrmat%AL_lagrange)
379 if(associated(fstrmat%AU_lagrange)) deallocate(fstrmat%AU_lagrange)
380 if(associated(fstrmat%Lagrange)) deallocate(fstrmat%Lagrange)
381
382 if( fstr_is_contact_active() ) then
383 allocate(fstrmat%AL_lagrange(ndof*fstrmat%numL_lagrange), stat=ierr)
384 if ( ierr /= 0 ) stop " Allocation error, fstrMAT%AL_lagrange "
385 fstrmat%AL_lagrange = 0.0d0
386 allocate(fstrmat%AU_lagrange(ndof*fstrmat%numU_lagrange), stat=ierr)
387 if ( ierr /= 0 ) stop " Allocation error, fstrMAT%AU_lagrange "
388 fstrmat%AU_lagrange = 0.0d0
389 allocate(fstrmat%Lagrange(fstrmat%num_lagrange))
390 fstrmat%Lagrange = 0.0d0
391 endif
392
393 if(associated(hecmat%B)) deallocate(hecmat%B)
394 allocate(hecmat%B(hecmat%NP*ndof+fstrmat%num_lagrange))
395 hecmat%B = 0.0d0
396
397 if(associated(hecmat%X)) deallocate(hecmat%X)
398 allocate(hecmat%X(hecmat%NP*ndof+fstrmat%num_lagrange))
399 hecmat%X = 0.0d0
400
401 if(associated(hecmat%D)) deallocate(hecmat%D)
402 allocate(hecmat%D(hecmat%NP*ndof**2+fstrmat%num_lagrange))
403 hecmat%D = 0.0d0
404 !
405 ! ---- For Parallel Contact with Multi-Partition Domains
406 if(paracontactflag.and.present(conmat)) then
407 if(associated(conmat%AL)) deallocate(conmat%AL)
408 allocate(conmat%AL(nn*conmat%NPL), stat=ierr)
409 if ( ierr /= 0 ) stop " Allocation error, conMAT%AL "
410 conmat%AL = 0.0d0
411
412 if(associated(conmat%AU)) deallocate(conmat%AU)
413 allocate(conmat%AU(nn*conmat%NPU), stat=ierr)
414 if ( ierr /= 0 ) stop " Allocation error, conMAT%AU "
415 conmat%AU = 0.0d0
416
417 if(associated(conmat%B)) deallocate(conmat%B)
418 allocate(conmat%B(conmat%NP*ndof+fstrmat%num_lagrange))
419 conmat%B = 0.0d0
420
421 if(associated(conmat%X)) deallocate(conmat%X)
422 allocate(conmat%X(conmat%NP*ndof+fstrmat%num_lagrange))
423 conmat%X = 0.0d0
424
425 if(associated(conmat%D)) deallocate(conmat%D)
426 allocate(conmat%D(conmat%NP*ndof**2+fstrmat%num_lagrange))
427 conmat%D = 0.0d0
428 endif
429
430 end subroutine constructnewmatrixstructure
431
432
434 subroutine insert_lagrange(i,id_lagrange,list_node,countNon0_lagrange)
435
436 type(noderelated) :: list_node
437 integer(kind=kint) :: i, id_lagrange
439 integer(kind=kint) :: countNon0_lagrange
441 integer(kind=kint) :: ierr, num_lagrange, location
442 integer(kind=kint) :: id_lagrange_save(1000)
443
444 character(len=1) :: answer
445
446 ierr = 0
447
448 num_lagrange = count(list_node%id_lagrange /= 0 )
449
450 ! if( i == 1 .and. num_lagrange /= 0) return
451 if( i == 1 .and. num_lagrange /= 0 .and. .not. permission) then
452 1 write(*,*) '##Error: node is both slave and master node simultaneously !'
453 write(*,*) ' Please check contact surface definition !'
454 write(*,'('' Do you want to continue(y/n)):'',$)')
455 read(*,'(A1)',err=1) answer
456 if(answer == 'Y' .OR. answer == 'y')then
457 permission = .true.
458 else
459 stop
460 endif
461 endif
462
463 if (num_lagrange == 0)then
464 list_node%num_lagrange = 1
465 list_node%id_lagrange(1) = id_lagrange
466 countnon0_lagrange = countnon0_lagrange + 1
467 else
468 id_lagrange_save(1:num_lagrange) = list_node%id_lagrange(1:num_lagrange)
469 location = find_locationinarray(id_lagrange,num_lagrange,list_node%id_lagrange)
470 if(location /= 0)then
471 num_lagrange = num_lagrange + 1
472 if( num_lagrange > size(list_node%id_lagrange)) then
473 deallocate(list_node%id_lagrange)
474 allocate(list_node%id_lagrange(num_lagrange),stat=ierr)
475 if( ierr /= 0 ) stop " Allocation error, list_nodeRelated%id_lagrange "
476 endif
477 list_node%num_lagrange = num_lagrange
478 list_node%id_lagrange(location) = id_lagrange
479 if(location /= 1) list_node%id_lagrange(1:location-1) = id_lagrange_save(1:location-1)
480 if(location /= num_lagrange) list_node%id_lagrange(location+1:num_lagrange) = id_lagrange_save(location:num_lagrange-1)
481 countnon0_lagrange = countnon0_lagrange + 1
482 endif
483 endif
484
485 end subroutine insert_lagrange
486
488 subroutine insert_node(id_node,list_node,countNon0_node)
489
490 type(noderelated) :: list_node
491 integer(kind=kint) :: id_node
493 integer(kind=kint) :: countNon0_node
494 integer(kind=kint) :: ierr, num_node, location
495 integer(kind=kint) :: id_node_save(1000)
496
497 ierr = 0
498
499 num_node = list_node%num_node
500
501 id_node_save(1:num_node) = list_node%id_node(1:num_node)
502 location = find_locationinarray(id_node,num_node,list_node%id_node)
503 if(location /= 0)then
504 num_node = num_node + 1
505 if( num_node > size(list_node%id_node)) then
506 deallocate(list_node%id_node)
507 allocate(list_node%id_node(num_node),stat=ierr)
508 if( ierr /= 0) stop " Allocation error, list_nodeRelated%id_node "
509 endif
510 list_node%num_node = num_node
511 list_node%id_node(location) = id_node
512 if(location /= 1) list_node%id_node(1:location-1) = id_node_save(1:location-1)
513 if(location /= num_node) list_node%id_node(location+1:num_node) = id_node_save(location:num_node-1)
514 countnon0_node = countnon0_node + 1
515 endif
516
517 end subroutine insert_node
518
519
521 integer function find_locationinarray(item,n,a)
522
523 integer(kind=kint) :: item, n
524 integer(kind=kint), pointer :: a(:)
525 integer(kind=kint) :: l, r, m
526
527 find_locationinarray = 0
528
529 l = 1 ; r = n ; m = (l+r)/2
530 if( item == a(l) .or. item == a(r) )then
531 return
532 elseif( item < a(l) )then
533 find_locationinarray = 1
534 return
535 elseif( item > a(r) )then
536 find_locationinarray = n + 1
537 return
538 endif
539
540 do while ( l <= r)
541 if( item > a(m) ) then
542 l = m + 1
543 m = (l + r)/2
544 elseif( item < a(m) ) then
545 r = m - 1
546 m = (l + r)/2
547 elseif( item == a(m) )then
548 return
549 endif
550 enddo
551
552 find_locationinarray = m + 1
553
554 end function find_locationinarray
555
556
558 subroutine reallocate_memory(num,list_node)
559
560 type(noderelated) :: list_node
561 integer(kind=kint) :: num
562 integer(kind=kint) :: num_node_org
564 integer(kind=kint) :: id_save(1000)
565 integer(kind=kint) :: ierr
566
567 num_node_org = size(list_node%id_node)
568 id_save(1:num_node_org) = list_node%id_node(1:num_node_org)
569 deallocate(list_node%id_node)
570 allocate(list_node%id_node(num_node_org+num),stat=ierr)
571 if( ierr /= 0) stop " reAllocation error, list_nodeRelated%id_node "
572 list_node%id_node = 0
573 list_node%id_node(1:num_node_org) = id_save(1:num_node_org)
574
575 end subroutine reallocate_memory
576
577
579 subroutine fstr_copy_lagrange_contact(fstrSOLID,fstrMAT)
580
581 type(fstr_solid) :: fstrSOLID
582 type(fstrst_matrix_contact_lagrange) :: fstrMAT
583 integer (kind=kint) :: id_lagrange, i, j
584
585 id_lagrange = 0
586
587 do i = 1, size(fstrsolid%contacts)
588 do j = 1, size(fstrsolid%contacts(i)%slave)
589 if( fstrsolid%contacts(i)%states(j)%state == contactfree ) cycle
590 id_lagrange = id_lagrange + 1
591 fstrmat%Lagrange(id_lagrange)=fstrsolid%contacts(i)%states(j)%multiplier(1)
592 enddo
593 enddo
594
595 end subroutine fstr_copy_lagrange_contact
596
597
599 logical function fstr_is_matrixstruct_symmetric(fstrSOLID,hecMESH)
600
601 type(fstr_solid ) :: fstrsolid
602 type(hecmwst_local_mesh) :: hecmesh
603 integer (kind=kint) :: is_in_contact
604
605 is_in_contact = 0
606 if( any(fstrsolid%contacts(:)%fcoeff /= 0.0d0) ) &
607 is_in_contact = 1
608 call hecmw_allreduce_i1(hecmesh, is_in_contact, hecmw_max)
609 fstr_is_matrixstruct_symmetric = (is_in_contact == 0)
610
612
613
This module encapsulate the basic functions of all elements provide by this software.
Definition: element.f90:43
integer, parameter fe_tri3n
Definition: element.f90:69
integer, parameter fe_quad4n
Definition: element.f90:72
This module provides functions of reconstructing.
subroutine getoriginallistofrelatednodes(np, num_lagrange)
Get original list of related nodes.
integer(kind=kint), save npl_org
subroutine constructnewmatrixstructure(hecmat, fstrmat, numnon0_node, numnon0_lagrange, conmat)
Construct new stiffness matrix structure.
integer(kind=kint), save npu_org
original number of non-zero items
subroutine fstr_save_originalmatrixstructure(hecmat)
This subroutine saves original matrix structure constructed originally by hecMW_matrix.
logical function fstr_is_matrixstruct_symmetric(fstrsolid, hecmesh)
this function judges whether sitiffness matrix is symmetric or not
subroutine fstr_mat_con_contact(cstep, hecmat, fstrsolid, fstrmat, infoctchange, conmat)
this subroutine reconstructs node-based (stiffness) matrix structure \corresponding to contact state
type(noderelated), dimension(:), pointer list_noderelated
current structure of matrix
subroutine fstr_copy_lagrange_contact(fstrsolid, fstrmat)
Copy Lagrange multipliers.
type(noderelated), dimension(:), pointer, save list_noderelated_org
original structure of matrix
subroutine getnewlistofrelatednodesandlagrangemultipliers(cstep, np, fstrsolid, countnon0lu_node, countnon0lu_lagrange)
Construct new list of related nodes and Lagrange multipliers. Here, a procedure similar to HEC_MW is ...
This module defined coomon data and basic structures for analysis.
Definition: m_fstr.f90:15
logical function fstr_iscontactactive(fstrsolid, nbc, cstep)
Definition: m_fstr.f90:1021
logical paracontactflag
PARALLEL CONTACT FLAG.
Definition: m_fstr.f90:84
This module provides functions to calcualte contact stiff matrix.
Definition: fstr_contact.f90:6
logical function fstr_is_contact_active()
Structure for Lagrange multiplier-related part of stiffness matrix (Lagrange multiplier-related matri...
Structure for defining stiffness matrix structure.