Rocstar  1.0
Rocstar multiphysics simulation application
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
MP/Source/UpdateStructuralSoln.f90
Go to the documentation of this file.
1 !*********************************************************************
2 !* Illinois Open Source License *
3 !* *
4 !* University of Illinois/NCSA *
5 !* Open Source License *
6 !* *
7 !* Copyright@2008, University of Illinois. All rights reserved. *
8 !* *
9 !* Developed by: *
10 !* *
11 !* Center for Simulation of Advanced Rockets *
12 !* *
13 !* University of Illinois *
14 !* *
15 !* www.csar.uiuc.edu *
16 !* *
17 !* Permission is hereby granted, free of charge, to any person *
18 !* obtaining a copy of this software and associated documentation *
19 !* files (the "Software"), to deal with the Software without *
20 !* restriction, including without limitation the rights to use, *
21 !* copy, modify, merge, publish, distribute, sublicense, and/or *
22 !* sell copies of the Software, and to permit persons to whom the *
23 !* Software is furnished to do so, subject to the following *
24 !* conditions: *
25 !* *
26 !* *
27 !* @ Redistributions of source code must retain the above copyright *
28 !* notice, this list of conditions and the following disclaimers. *
29 !* *
30 !* @ Redistributions in binary form must reproduce the above *
31 !* copyright notice, this list of conditions and the following *
32 !* disclaimers in the documentation and/or other materials *
33 !* provided with the distribution. *
34 !* *
35 !* @ Neither the names of the Center for Simulation of Advanced *
36 !* Rockets, the University of Illinois, nor the names of its *
37 !* contributors may be used to endorse or promote products derived *
38 !* from this Software without specific prior written permission. *
39 !* *
40 !* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, *
41 !* EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES *
42 !* OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND *
43 !* NONINFRINGEMENT. IN NO EVENT SHALL THE CONTRIBUTORS OR *
44 !* COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER *
45 !* LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, *
46 !* ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE *
47 !* USE OR OTHER DEALINGS WITH THE SOFTWARE. *
48 !*********************************************************************
49 !* Please acknowledge The University of Illinois Center for *
50 !* Simulation of Advanced Rockets in works and publications *
51 !* resulting from this software or its derivatives. *
52 !*********************************************************************
53 
55 
56  USE rocstar_rocfrac
58 
59 CONTAINS
60 
61  SUBROUTINE updatestructural(glb,NumProcs,Rnet)
62 
63  IMPLICIT NONE
64 
65  include 'mpif.h'
66 
67  TYPE(rocfrac_global) :: glb
68 
69  INTEGER :: numnp ! Number of Node Points
70  INTEGER :: numelvol ! Number of Volumetric Elements
71  INTEGER :: nummatvol ! Number of Volumetreic Materials
72  INTEGER :: nummatcoh ! Number of Cohesive Materials
73  INTEGER :: numprocs ! Number of processors
74  INTEGER :: ieltype ! Order of element (4:4node, 5:4nodeEnhanced, 10:10node)
75  REAL*8, DIMENSION(1:3*glb%NumNP) :: rnet
76  REAL*8, POINTER, DIMENSION(:,:) :: meshcoor
77  INTEGER, POINTER, DIMENSION(:) :: matidvol
78  INTEGER, POINTER, DIMENSION(:,:) :: elconnvol
79  REAL*8, POINTER, DIMENSION(:) :: e, xnu, rho
80  REAL*8, POINTER, DIMENSION(:) :: disp ! Nodal Displacement
81  REAL*8, POINTER, DIMENSION(:) :: deltan, deltat
82  REAL*8, POINTER, DIMENSION(:) :: velohalf ! velocity
83  REAL*8, POINTER, DIMENSION(:) :: velobar, accelbar
84  REAL*8 :: kappadamp
85  LOGICAL :: aleenabled, dampenabled
86  INTEGER, POINTER, DIMENSION(:) :: isolntype
87 ! LOGICAL :: DefConfig
88  INTEGER :: totnumneighprocs
89  INTEGER :: mpi_comm_rocfrac
90 !-- Non-block receive, Non-block send request arrays
91  INTEGER, POINTER, DIMENSION(:) :: reqrcv, reqsnd
92  INTEGER, POINTER, DIMENSION(:,:) :: statrcv, statsnd
93  INTEGER :: totnumndcomm
94 
95  INTEGER :: numelpartbndry
96  INTEGER, DIMENSION(:), POINTER :: numelpartbndrymat
97  INTEGER, DIMENSION(:), POINTER :: neighproclist
98  INTEGER, DIMENSION(:), POINTER :: numndcomm
99  INTEGER, DIMENSION(:), POINTER :: numelvolmat
100 
101  TYPE(send_buf), DIMENSION(:), POINTER :: ndcommlist
102 
103  TYPE(rcv_buf), DIMENSION(:), POINTER :: recvdatafrm
104 
105  TYPE(send_buf), POINTER :: pndcomm
106  TYPE(rcv_buf), POINTER :: precvdf
107 
108  REAL*8, ALLOCATABLE, DIMENSION(:) :: buf
109 
110  INTEGER :: j, j1, k , k1, k2
111  INTEGER :: ierr,ianalysistype
112  INTEGER :: elemstart, elemend
113  REAL*8, POINTER, DIMENSION(:) :: sigmamax, taumax, sinit
114  REAL*8, POINTER, DIMENSION(:,:) :: sthresh1 ,sthresh2
115 
116 !------ CALCULATE THE SUBMESH'S BOUNDARY INTERNAL FORCE VECTOR
117 
118 
119  numnp = glb%NumNP
120  numelvol = glb%NumElVol
121  nummatvol = glb%NumMatVol
122  nummatcoh = glb%NumMatCoh
123  ieltype = glb%iEltype
124 
125 ! Rnet => glb%Rnet
126  meshcoor => glb%MeshCoor
127  matidvol => glb%MatIdVol
128  elconnvol => glb%ElConnVol
129  disp => glb%Disp
130  deltan => glb%deltan
131  deltat => glb%deltat
132  sigmamax => glb%SigmaMax
133  taumax => glb%TauMax
134  sthresh1 => glb%Sthresh1
135  sthresh2 => glb%Sthresh2
136  sinit => glb%Sinit
137  velohalf => glb%VeloHalf
138  velobar => glb%VeloBar
139  accelbar => glb%AccelBar
140  kappadamp = glb%KappaDamp
141  aleenabled = glb%ALEenabled
142  dampenabled = glb%DampEnabled
143  isolntype => glb%iSolnType
144 ! DefConfig = glb%DefConfig
145  totnumneighprocs = glb%TotNumNeighProcs
146  mpi_comm_rocfrac = glb%MPI_COMM_ROCFRAC
147  reqrcv => glb%ReqRcv
148  reqsnd => glb%ReqSnd
149  statrcv => glb%StatRcv
150  statsnd => glb%StatSnd
151  totnumndcomm = glb%TotNumNdComm
152  numelpartbndry = glb%NumElPartBndry
153  numelpartbndrymat => glb%NumElPartBndryMat
154  neighproclist => glb%NeighProcList
155  numndcomm => glb%NumNdComm
156  numelvolmat => glb%NumElVolMat
157  ndcommlist => glb%NdCommList
158  recvdatafrm => glb%RecvDataFrm
159  e => glb%E
160  xnu => glb%xnu
161  rho => glb%rho
162 
163  elemstart = 1
164 
165  DO j = 1, nummatvol
166 
167  elemend = numelpartbndrymat(j) + elemstart - 1
168 
169  ianalysistype = isolntype(j)
170 
171  IF(ieltype.EQ.4)THEN
172  CALL internalforce_v3d4( glb, rnet, elemstart, elemend, ianalysistype)
173  ELSE IF(ieltype.EQ.10)THEN
174  CALL internalforce_v3d10( glb, rnet, elemstart, elemend, ianalysistype)
175  ELSE IF(ieltype.EQ.8)THEN
176  CALL internalforce_v3d8( glb, rnet, elemstart, elemend, ianalysistype)
177  ENDIF
178 
179  elemstart = numelvolmat(j) + elemstart
180  ENDDO
181 
182  DO j = 1, nummatcoh
183 
184  ! Transfer the cohesive solution to the 'Top' Surface
185 
186 ! PRINT*,'finished readsdv23'
187 ! CALL MPI_BARRIER(glb%MPI_COMM_ROCFRAC,ierr)
188 ! STOP
189 
190 
191  CALL c3d6nm( glb%nsubn1, glb%nsubf1, &
192  glb%nsubn2, glb%nsubf2, &
193  glb%sd_coor1, glb%sd_subfaces1, &
194  glb%sd_subface_parents1, &
195  glb%sd_subface_parents2, &
196  glb%sd_subface_counterparts1, &
197  glb%sd_subface_nat_coors1, &
198  glb%sd_subface_nat_coors2, &
199  glb%FaceOfVolEl1, &
200  glb%FaceOfVolEl2, &
201  glb%NumNp, glb%NumElVol, &
202  glb%ElConnVol, &
203  glb%nf1,glb%nf2, &
204  glb%MapFaceEl2Vol1, &
205  glb%MapFaceEl2Vol2, &
206  deltan, deltat, sigmamax, taumax, sinit, &
207  rnet, disp, sthresh1, glb%NumMatCoh,j,-1.d0)
208 
209  ! Transfer the cohesive solution to the 'Bottom' Surface
210 
211 !!$ CALL c3d6nm( glb%nsubn2, glb%nsubf2, &
212 !!$ glb%sd_coor2, glb%sd_subfaces2, &
213 !!$ glb%sd_subface_parents2, &
214 !!$ glb%sd_subface_counterparts2, &
215 !!$ glb%sd_subface_nat_coors2, &
216 !!$ glb%FaceOfVolEl2, &
217 !!$ glb%NumNp, glb%NumElVol, &
218 !!$ glb%ElConnVol, &
219 !!$ glb%nf2, &
220 !!$ glb%MapFaceEl2Vol2, &
221 !!$ deltan, deltat, SigmaMax, TauMax, Sinit, &
222 !!$ Rnet, Disp, Sthresh1, glb%NumMatCoh,j,1.d0)
223 
224 
225  ENDDO
226 
227 !
228 !----- FORM THE BUFFER CONTAINING COMMUNICATED NODAL VALUES
229 !
230  ALLOCATE(buf(1:totnumndcomm))
231  k1 = 1
232  DO j1 = 1, totnumneighprocs
233  k = neighproclist(j1)
234  pndcomm => ndcommlist(j1)
235  DO j = 1, numndcomm(j1)
236  k2 = 3*pndcomm%NdId(j)
237  buf(k1) = rnet( k2 - 2 )
238  buf(k1+1) = rnet( k2 - 1 )
239  buf(k1+2) = rnet( k2 )
240  k1 = k1 + 3
241  ENDDO
242  ENDDO
243 
244 !
245 !-MPI- RECEIVE THE INTERNAL FORCE VECTOR FROM THE NEIGHBORS
246 !
247  DO j1 = 1, totnumneighprocs
248  k = neighproclist(j1)+1
249  CALL mpi_irecv(recvdatafrm(k)%rcvbuf(1),numndcomm(j1)*3, &
250  mpi_double_precision,k-1,10,mpi_comm_rocfrac, &
251  reqrcv(j1),ierr)
252  ENDDO
253 !
254 !-MPI- SEND THE INTERNAL FORCE VECTOR TO THE NEIGHBORS
255 !
256  k2 = 1
257  DO j1 = 1, totnumneighprocs
258  k = neighproclist(j1)
259  CALL mpi_isend(buf(k2), numndcomm(j1)*3,mpi_double_precision, &
260  k,10,mpi_comm_rocfrac,reqsnd(j1),ierr)
261  k2 = k2 + numndcomm(j1)*3
262  ENDDO
263 
264 !
265 !-MPI- WAIT FOR INTERNAL FORCE VECTOR COMMUNICATION TO COMPLETE
266 !
267 
268  IF(aleenabled)THEN
269  IF(ieltype.EQ.4)THEN
270  CALL v3d4_ale(velobar,accelbar,disp,velohalf,rnet, &
271  e,xnu,rho,numnp,nummatvol, &
272  numelvol,matidvol,elconnvol,meshcoor, &
273  numelpartbndry+1,numelvol)
274  ELSE
275  CALL v3d10_ale(velobar,accelbar,disp,velohalf,rnet, &
276  e,xnu,rho,numnp,nummatvol, &
277  numelvol,matidvol,elconnvol,meshcoor, &
278  numelpartbndry+1,numelvol)
279  ENDIF
280  ENDIF
281 
282 !-- (11) calculate R_in, R_damp
283 
284 !------ CALCULATE THE SUBMESH'S BOUNDARY INTERNAL FORCE VECTOR
285 
286  elemend = 0
287  DO j = 1, nummatvol
288 
289  elemstart = elemend + numelpartbndrymat(j) + 1
290 
291  elemend = numelvolmat(j) + elemend
292 
293  IF(ieltype.EQ.4)THEN
294  CALL internalforce_v3d4( glb, rnet, elemstart, elemend, ianalysistype)
295  ELSE IF(ieltype.EQ.10)THEN
296  CALL internalforce_v3d10( glb, rnet, elemstart, elemend, ianalysistype)
297  ELSE IF(ieltype.EQ.8)THEN
298  CALL internalforce_v3d8( glb, rnet, elemstart, elemend, ianalysistype)
299  ENDIF
300  ENDDO
301 
302  IF(totnumneighprocs.GT.0)THEN
303  CALL mpi_waitall(totnumneighprocs,reqrcv,statrcv,ierr)
304  CALL mpi_waitall(totnumneighprocs,reqsnd,statsnd,ierr)
305  ENDIF
306 
307  DEALLOCATE(buf)
308 
309 !
310 !----- ADD NEIGHBOR'S CONTRIBUTION TO THE INTERNAL FORCE VECTOR
311 !
312  DO j1 = 1, totnumneighprocs
313  k = neighproclist(j1)+1
314  k1 = 1
315  pndcomm => ndcommlist(j1)
316  precvdf => recvdatafrm(k)
317  DO j = 1, numndcomm(j1)
318  k2 = ( pndcomm%NdId(j) )*3
319  rnet(k2-2)= rnet(k2-2) + precvdf%rcvbuf(k1)
320  rnet(k2-1)= rnet(k2-1) + precvdf%rcvbuf(k1+1)
321  rnet(k2) = rnet(k2) + precvdf%rcvbuf(k1+2)
322  k1 = k1 + 3
323  ENDDO
324  ENDDO
325 
326  RETURN
327  END SUBROUTINE updatestructural
328 
329  SUBROUTINE updatestructuralht(glb,NumProcs,Rnet,RnetHT)
330 
331  IMPLICIT NONE
332 
333  include 'mpif.h'
334 
335  TYPE(rocfrac_global) :: glb
336 
337  INTEGER :: numnp ! Number of Node Points
338  INTEGER :: numelvol ! Number of Volumetric Elements
339  INTEGER :: nummatvol ! Number of Volumetreic Materials
340  INTEGER :: numprocs ! Number of processors
341  INTEGER :: ieltype ! Order of element (4:4node, 5:4nodeEnhanced, 10:10node)
342  REAL*8, DIMENSION(1:3*glb%NumNP) :: rnet
343  REAL*8, DIMENSION(1:glb%NumNP) :: rnetht
344  REAL*8, POINTER, DIMENSION(:,:) :: meshcoor
345  INTEGER, POINTER, DIMENSION(:) :: matidvol
346  INTEGER, POINTER, DIMENSION(:,:) :: elconnvol
347  REAL*8, POINTER, DIMENSION(:) :: e, xnu, rho, cp, kappaht
348  REAL*8, POINTER, DIMENSION(:) :: disp ! Nodal Displacement
349  REAL*8, POINTER, DIMENSION(:) :: velohalf ! velocity
350  REAL*8, POINTER, DIMENSION(:) :: velobar, accelbar
351  REAL*8 :: kappadamp
352  LOGICAL :: aleenabled, dampenabled
353  INTEGER, POINTER, DIMENSION(:) :: isolntype
354 ! LOGICAL :: DefConfig
355  INTEGER :: totnumneighprocs
356  INTEGER :: mpi_comm_rocfrac
357 !-- Non-block receive, Non-block send request arrays
358  INTEGER, POINTER, DIMENSION(:) :: reqrcv, reqsnd
359  INTEGER, POINTER, DIMENSION(:,:) :: statrcv, statsnd
360  INTEGER :: totnumndcomm
361 
362  INTEGER :: numelpartbndry
363  INTEGER, DIMENSION(:), POINTER :: numelpartbndrymat
364  INTEGER, DIMENSION(:), POINTER :: neighproclist
365  INTEGER, DIMENSION(:), POINTER :: numndcomm
366  INTEGER, DIMENSION(:), POINTER :: numelvolmat
367 
368  integer :: k4
369 
370  TYPE(send_buf), DIMENSION(:), POINTER :: ndcommlist
371 
372  TYPE(rcv_buf), DIMENSION(:), POINTER :: recvdatafrm
373 
374  TYPE(send_buf), POINTER :: pndcomm
375  TYPE(rcv_buf), POINTER :: precvdf
376 
377  REAL*8, ALLOCATABLE, DIMENSION(:) :: buf
378 
379  INTEGER :: i, j, j1, k , k1, k2
380  INTEGER :: ierr,ianalysistype
381  INTEGER :: elemstart, elemend
382  logical :: heattranssoln
383 
384  integer :: numndsbcht
385  real*8, DIMENSION(:), POINTER :: temperature
386  integer, DIMENSION(:,:), POINTER :: bcflaght
387  real*8, DIMENSION(:,:), POINTER :: bcvalueht
388 
389 !------ CALCULATE THE SUBMESH'S BOUNDARY INTERNAL FORCE VECTOR
390 
391 
392  numnp = glb%NumNP
393  numelvol = glb%NumElVol
394  nummatvol = glb%NumMatVol
395  ieltype = glb%iEltype
396 
397 ! Rnet => glb%Rnet
398  meshcoor => glb%MeshCoor
399  matidvol => glb%MatIdVol
400  elconnvol => glb%ElConnVol
401  disp => glb%Disp
402  velohalf => glb%VeloHalf
403  velobar => glb%VeloBar
404  accelbar => glb%AccelBar
405  kappadamp = glb%KappaDamp
406  aleenabled = glb%ALEenabled
407  dampenabled = glb%DampEnabled
408  isolntype => glb%iSolnType
409 ! DefConfig = glb%DefConfig
410  totnumneighprocs = glb%TotNumNeighProcs
411  mpi_comm_rocfrac = glb%MPI_COMM_ROCFRAC
412  reqrcv => glb%ReqRcv
413  reqsnd => glb%ReqSnd
414  statrcv => glb%StatRcv
415  statsnd => glb%StatSnd
416  totnumndcomm = glb%TotNumNdComm
417  numelpartbndry = glb%NumElPartBndry
418  numelpartbndrymat => glb%NumElPartBndryMat
419  neighproclist => glb%NeighProcList
420  numndcomm => glb%NumNdComm
421  numelvolmat => glb%NumElVolMat
422  ndcommlist => glb%NdCommList
423  recvdatafrm => glb%RecvDataFrm
424  e => glb%E
425  xnu => glb%xnu
426  rho => glb%rho
427  cp => glb%Cp
428  kappaht => glb%KappaHT
429  heattranssoln = glb%HeatTransSoln
430 
431  elemstart = 1
432 
433  numndsbcht = glb%NumNdsBCHT
434  bcflaght => glb%BCFlagHT
435  temperature => glb%Temperature
436  bcvalueht => glb%BCvalueHT
437  DO j = 1, nummatvol
438 
439  elemend = numelpartbndrymat(j) + elemstart - 1
440 
441  ianalysistype = isolntype(j)
442 
443 
444  DO i = 1,glb%NumNdsBCHT
445  k4 = glb%BCFlagHT(1,i) ! node
446 
447  IF (glb%BCFlagHT(2,i).EQ.0) THEN ! impose temperature
448  glb%Temperature(k4) = glb%BCvalueHT(1,i) ! *glb%prop
449 
450  ENDIF
451  ENDDO
452 
453  IF(ieltype.EQ.4)THEN
454 
455  glb%S11(:,:) = 0.d0
456  glb%S22(:,:) = 0.d0
457  glb%S33(:,:) = 0.d0
458 
459  CALL v3d4_thermal(numelvol, numnp, elconnvol, meshcoor, kappaht, &
460  rnetht, temperature, rho, cp, matidvol, nummatvol,velobar,elemstart, elemend)
461 
462 
463  glb%Temperature(:) = glb%Temperature(:) + glb%DT*glb%CapctInv(:)*rnetht(:)
464 
465  DO i = 1,glb%NumNdsBCHT
466  k4 = glb%BCFlagHT(1,i) ! node
467 
468  IF (glb%BCFlagHT(2,i).EQ.0) THEN ! impose temperature
469  glb%Temperature(k4) = glb%BCvalueHT(1,i) ! *glb%prop
470  ENDIF
471  ENDDO
472 
473  CALL internalforce_v3d4ht( glb, rnet, elemstart, elemend, ianalysistype, temperature)
474 
475 !!$ CALL v3d4_thermalExp(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,Rnet,glb%ci, &
476 !!$ glb%S11,glb%S22,glb%S33, &
477 !!$ glb%NumNP,ElemStart,ElemEnd,glb%NumElVol,glb%NumMatVol, glb%alpha, Temperature, glb%Temperature0)
478 
479  ELSE
480 
481  glb%S11(:,:) = 0.d0
482  glb%S22(:,:) = 0.d0
483  glb%S33(:,:) = 0.d0
484 
485  CALL internalforce_v3d10( glb, rnet, elemstart, elemend, ianalysistype)
486 
487  Call v3d10_thermal(numelvol, numnp, elconnvol, meshcoor, kappaht, &
488  rnetht, temperature, rho, cp, matidvol, nummatvol, velobar,elemstart, elemend)
489 
490  glb%Temperature(:) = glb%Temperature(:) + glb%DT*glb%CapctInv(:)*rnetht(:)
491 
492  DO i = 1,glb%NumNdsBCHT
493  k4 = glb%BCFlagHT(1,i) ! node
494 
495  IF (glb%BCFlagHT(2,i).EQ.0) THEN ! impose temperature
496  glb%Temperature(k4) = glb%BCvalueHT(1,i) !*glb%prop
497  ENDIF
498  ENDDO
499 
500 
501  CALL v3d10_thermalexp(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,rnet,glb%ci, &
502  glb%S11,glb%S22,glb%S33, &
503  glb%NumNP,elemstart,elemend,glb%NumElVol,glb%NumMatVol, glb%alpha, temperature, glb%Temperature0)
504 
505  ENDIF
506 
507  elemstart = numelvolmat(j) + elemstart
508  ENDDO
509 !
510 !----- FORM THE BUFFER CONTAINING COMMUNICATED NODAL VALUES
511 !
512  ALLOCATE(buf(1:totnumndcomm+totnumndcomm/3))
513  k1 = 1
514  DO j1 = 1, totnumneighprocs
515  k = neighproclist(j1)
516  pndcomm => ndcommlist(j1)
517  DO j = 1, numndcomm(j1)
518  k2 = 3*pndcomm%NdId(j)
519  buf(k1) = rnet( k2 - 2 )
520  buf(k1+1) = rnet( k2 - 1 )
521  buf(k1+2) = rnet( k2 )
522  buf(k1+3) = rnetht(pndcomm%NdId(j))
523  k1 = k1 + 4
524  ENDDO
525  ENDDO
526 
527 !
528 !-MPI- RECEIVE THE INTERNAL FORCE VECTOR FROM THE NEIGHBORS
529 !
530  DO j1 = 1, totnumneighprocs
531  k = neighproclist(j1)+1
532  CALL mpi_irecv(recvdatafrm(k)%rcvbuf(1),numndcomm(j1)*4, &
533  mpi_double_precision,k-1,10,mpi_comm_rocfrac, &
534  reqrcv(j1),ierr)
535  ENDDO
536 !
537 !-MPI- SEND THE INTERNAL FORCE VECTOR TO THE NEIGHBORS
538 !
539  k2 = 1
540  DO j1 = 1, totnumneighprocs
541  k = neighproclist(j1)
542  CALL mpi_isend(buf(k2), numndcomm(j1)*4,mpi_double_precision, &
543  k,10,mpi_comm_rocfrac,reqsnd(j1),ierr)
544  k2 = k2 + numndcomm(j1)*4
545  ENDDO
546 
547 !
548 !-MPI- WAIT FOR INTERNAL FORCE VECTOR COMMUNICATION TO COMPLETE
549 !
550 
551  IF(aleenabled)THEN
552  IF(ieltype.EQ.4)THEN
553  CALL v3d4_ale(velobar,accelbar,disp,velohalf,rnet, &
554  e,xnu,rho,numnp,nummatvol, &
555  numelvol,matidvol,elconnvol,meshcoor, &
556  numelpartbndry+1,numelvol)
557  ELSE
558  CALL v3d10_ale(velobar,accelbar,disp,velohalf,rnet, &
559  e,xnu,rho,numnp,nummatvol, &
560  numelvol,matidvol,elconnvol,meshcoor, &
561  numelpartbndry+1,numelvol)
562  ENDIF
563  ENDIF
564 
565 !-- (11) calculate R_in, R_damp
566 
567 !------ CALCULATE THE SUBMESH'S BOUNDARY INTERNAL FORCE VECTOR
568 
569  elemend = 0
570  DO j = 1, nummatvol
571 
572  elemstart = elemend + numelpartbndrymat(j) + 1
573 
574  elemend = numelvolmat(j) + elemend
575 
576  DO i = 1,glb%NumNdsBCHT
577  k4 = glb%BCFlagHT(1,i) ! node
578 
579  IF (glb%BCFlagHT(2,i).EQ.0) THEN ! impose temperature
580  glb%Temperature(k4) = glb%BCvalueHT(1,i) !*glb%prop
581 
582  ENDIF
583  ENDDO
584 
585  IF(ieltype.EQ.4)THEN
586 
587  CALL v3d4_thermal(numelvol, numnp, elconnvol, meshcoor, kappaht, &
588  rnetht, temperature, rho, cp, matidvol, nummatvol,velobar,elemstart, elemend)
589 
590  glb%Temperature(:) = glb%Temperature(:) + glb%DT*glb%CapctInv(:)*rnetht(:)
591 
592  DO i = 1,glb%NumNdsBCHT
593  k4 = glb%BCFlagHT(1,i) ! node
594 
595  IF (glb%BCFlagHT(2,i).EQ.0) THEN ! impose temperature
596  glb%Temperature(k4) = glb%BCvalueHT(1,i) ! *glb%prop
597  ENDIF
598  ENDDO
599 
600  CALL internalforce_v3d4ht( glb, rnet, elemstart, elemend, ianalysistype, temperature)
601 
602 !!$ CALL v3d4_thermalExp(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,Rnet,glb%ci, &
603 !!$ glb%S11,glb%S22,glb%S33, &
604 !!$ glb%NumNP,ElemStart,ElemEnd,glb%NumElVol,glb%NumMatVol, glb%alpha, Temperature, glb%Temperature0)
605 
606  ELSE
607  CALL internalforce_v3d10( glb, rnet, elemstart, elemend, ianalysistype)
608 
609  Call v3d10_thermal(numelvol, numnp, elconnvol, meshcoor, kappaht, &
610  rnetht, temperature, rho, cp, matidvol, nummatvol, velobar,elemstart, elemend)
611 
612  glb%Temperature(:) = glb%Temperature(:) + glb%DT*glb%CapctInv(:)*rnetht(:)
613 
614  DO i = 1,glb%NumNdsBCHT
615  k4 = glb%BCFlagHT(1,i) ! node
616 
617  IF (glb%BCFlagHT(2,i).EQ.0) THEN ! impose temperature
618  glb%Temperature(k4) = glb%BCvalueHT(1,i) !*glb%prop
619  ENDIF
620  ENDDO
621 
622  CALL v3d10_thermalexp(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,rnet,glb%ci, &
623  glb%S11,glb%S22,glb%S33, &
624  glb%NumNP,elemstart,elemend,glb%NumElVol,glb%NumMatVol, glb%alpha, temperature, glb%Temperature0)
625  ENDIF
626  ENDDO
627 
628  IF(totnumneighprocs.GT.0)THEN
629  CALL mpi_waitall(totnumneighprocs,reqrcv,statrcv,ierr)
630  CALL mpi_waitall(totnumneighprocs,reqsnd,statsnd,ierr)
631  ENDIF
632 
633  DEALLOCATE(buf)
634 
635 !
636 !----- ADD NEIGHBOR'S CONTRIBUTION TO THE INTERNAL FORCE VECTOR
637 !
638  DO j1 = 1, totnumneighprocs
639  k = neighproclist(j1)+1
640  k1 = 1
641  pndcomm => ndcommlist(j1)
642  precvdf => recvdatafrm(k)
643  DO j = 1, numndcomm(j1)
644  k2 = ( pndcomm%NdId(j) )*3
645  rnet(k2-2)= rnet(k2-2) + precvdf%rcvbuf(k1)
646  rnet(k2-1)= rnet(k2-1) + precvdf%rcvbuf(k1+1)
647  rnet(k2) = rnet(k2) + precvdf%rcvbuf(k1+2)
648  rnetht(pndcomm%NdId(j)) = rnetht(pndcomm%NdId(j)) + precvdf%rcvbuf(k1+3)
649  k1 = k1 + 4
650  ENDDO
651  ENDDO
652 
653  RETURN
654  END SUBROUTINE updatestructuralht
655 
656  SUBROUTINE internalforce_v3d4( glb, Rnet, ElemStart, ElemEnd, iAnalysisType)
657 
658  IMPLICIT NONE
659 
660  TYPE(rocfrac_global) :: glb
661  REAL*8, DIMENSION(1:3*glb%NumNP) :: rnet
662 
663  INTEGER, INTENT(IN) :: elemstart, elemend, ianalysistype
664 
665  IF (ianalysistype.EQ.0) THEN
666  CALL v3d4_nl_arruda_boyce(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,rnet,glb%Disp,&
667  glb%S11,glb%S22,glb%S33,glb%S12,glb%S23,glb%S13, &
668  glb%NumNP,elemstart,elemend,glb%NumElVol,glb%NumMatVol, &
669  glb%xmu,glb%xkappa)
670  ELSE IF(ianalysistype.EQ.1)THEN
671  CALL v3d4_nl(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,rnet,glb%Disp,glb%ci, &
672  glb%S11,glb%S22,glb%S33,glb%S12,glb%S23,glb%S13, &
673  glb%NumNP,elemstart,elemend,glb%NumElVol,glb%NumMatVol,glb%xmu,glb%xlambda)
674  ELSE IF(ianalysistype.EQ.-1)THEN
675  CALL v3d4_neohookeanincompress(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,rnet,glb%Disp,&
676  glb%S11,glb%S22,glb%S33,glb%S12,glb%S23,glb%S13, &
677  glb%NumNP,elemstart,elemend,glb%NumElVol,glb%NumMatVol, &
678  glb%xmu,glb%xkappa)
679  ELSE
680 !------- Linear elastic
681  CALL v3d4(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,rnet,glb%Disp,glb%ci, &
682  glb%S11,glb%S22,glb%S33,glb%S12,glb%S23,glb%S13, &
683  glb%NumNP,elemstart,elemend,glb%NumElVol,glb%NumMatVol)
684  ENDIF
685 
686 ! Damping Rdamp
687 
688  IF(glb%DampEnabled)THEN
689  CALL v3d4_damping(glb%VeloHalf,rnet,glb%NumNP, &
690  glb%NumElVol,glb%ElConnVol,glb%MeshCoor, &
691  elemstart,elemend, glb%KappaDamp)
692  END IF
693 
694 
695  END SUBROUTINE internalforce_v3d4
696 
697  SUBROUTINE internalforce_v3d4ht( glb, Rnet, ElemStart, ElemEnd, iAnalysisType, Temperature)
698 
699  IMPLICIT NONE
700 
701  TYPE(rocfrac_global) :: glb
702  REAL*8, DIMENSION(1:3*glb%NumNP) :: rnet
703  REAL*8, DIMENSION(1:glb%NumNP) :: temperature
704 
705  INTEGER, INTENT(IN) :: elemstart, elemend, ianalysistype
706 
707  IF (ianalysistype.EQ.0) THEN
708  CALL v3d4_nl_arruda_boyce(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,rnet,glb%Disp,&
709  glb%S11,glb%S22,glb%S33,glb%S12,glb%S23,glb%S13, &
710  glb%NumNP,elemstart,elemend,glb%NumElVol,glb%NumMatVol, &
711  glb%xmu,glb%xkappa)
712  ELSE IF(ianalysistype.EQ.1)THEN
713  CALL v3d4_nl(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,rnet,glb%Disp,glb%ci, &
714  glb%S11,glb%S22,glb%S33,glb%S12,glb%S23,glb%S13, &
715  glb%NumNP,elemstart,elemend,glb%NumElVol,glb%NumMatVol,glb%xmu,glb%xlambda)
716  ELSE IF(ianalysistype.EQ.-1)THEN
717  CALL v3d4_neohookeanincompress(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,rnet,glb%Disp,&
718  glb%S11,glb%S22,glb%S33,glb%S12,glb%S23,glb%S13, &
719  glb%NumNP,elemstart,elemend,glb%NumElVol,glb%NumMatVol, &
720  glb%xmu,glb%xkappa)
721  ELSE
722 !------- Linear elastic
723  CALL v3d4_thermalexp2(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,rnet,glb%Disp,glb%ci, &
724  glb%S11,glb%S22,glb%S33,glb%S12,glb%S23,glb%S13, &
725  glb%NumNP,elemstart,elemend,glb%NumElVol,glb%NumMatVol, &
726  glb%alpha,temperature,glb%Temperature0)
727  ENDIF
728 
729 ! Damping Rdamp
730 
731  IF(glb%DampEnabled)THEN
732  CALL v3d4_damping(glb%VeloHalf,rnet,glb%NumNP, &
733  glb%NumElVol,glb%ElConnVol,glb%MeshCoor, &
734  elemstart,elemend, glb%KappaDamp)
735  END IF
736 
737 
738  END SUBROUTINE internalforce_v3d4ht
739 
740 
741  SUBROUTINE internalforce_v3d10( glb, Rnet, ElemStart, ElemEnd, iAnalysisType)
742 
743  IMPLICIT NONE
744 
745  TYPE(rocfrac_global) :: glb
746  REAL*8, DIMENSION(1:3*glb%NumNP) :: rnet
747 
748  INTEGER, INTENT(IN) :: elemstart, elemend, ianalysistype
749 
750  IF(glb%DebondPart)THEN
751 
752  CALL v3d10_nl_huang(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,rnet,glb%Disp, &
753  glb%S11,glb%S22,glb%S33,glb%S12,glb%S23,glb%S13, &
754  glb%NumNP,elemstart,elemend,glb%NumElVol,glb%NumMatVol, &
755  glb%STATEV_Part1,glb%STATEV_Part2,glb%NSTATEV,glb%MATRIX,glb%NMATRIX, &
756  glb%PARTICLE,glb%NPARTICLE,glb%NPARTICLETYPE,glb%INTERFAC,glb%NINTERFAC,glb%StrainTrace)
757 
758  ELSE IF(glb%DebondPart_Matous)THEN
759 
760 
761  IF(glb%Debug_State) print*,'Starting v3d10_nl_matous', elemstart,elemend
762 
763  CALL v3d10_nl_matous(glb%MeshCoor,glb%ElConnVol,rnet,glb%Disp, &
764  glb%S11,glb%S22,glb%S33,glb%S12,glb%S23,glb%S13, &
765  glb%NumNP,elemstart,elemend,glb%NumElVol, &
766  glb%p1, glb%p2, glb%Yin, glb%SoftParam, glb%BulkMod(2), glb%BulkMod(1), &
767  glb%ShrMod(2), glb%ShrMod(1), glb%PoisRat(2), glb%PoisRat(1), &
768  glb%a_eta, glb%a_zeta,&
769  glb%cm, glb%c2, glb%cd, glb%Lo, glb%L_tensor(:,:,1), &
770  glb%L_tensor(:,:,2), glb%M_tensor(:,:,1), glb%M_tensor(:,:,2), glb%L_bar, &
771  glb%StrainOld)
772 
773  IF(glb%Debug_State) print*,'Finished v3d10_nl_matous', elemstart,elemend
774 
775 
776  ELSE
777 
778 
779  IF (ianalysistype.EQ.0) THEN
780 
781  IF ( glb%ArtificialDamping)THEN
782 
783  CALL v3d10_nl_arruda_boyce_damping(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,rnet,glb%Disp, &
784  glb%S11,glb%S22,glb%S33,glb%S12,glb%S23,glb%S13, &
785  glb%NumNP,elemstart,elemend,glb%NumElVol,glb%NumMatVol, &
786  glb%xmu,glb%xkappa, glb%rho, glb%cd_fastest, glb%DetF_old, glb%VeloHalf, glb%Dt)
787 
788  ELSE
789 
790  CALL v3d10_nl_arruda_boyce(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,rnet,glb%Disp, &
791  glb%S11,glb%S22,glb%S33,glb%S12,glb%S23,glb%S13, &
792  glb%NumNP,elemstart,elemend,glb%NumElVol,glb%NumMatVol, &
793  glb%xmu,glb%xkappa)
794  ENDIF
795 
796  ELSE IF (ianalysistype.EQ.1)THEN
797 
798  IF(glb%iElIntgratn.EQ.0)THEN
799 
800  IF ( glb%ArtificialDamping)THEN
801 
802 
803  CALL v3d10_nl_damping(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,rnet,glb%Disp,glb%ci, &
804  glb%S11,glb%S22,glb%S33,glb%S12,glb%S23,glb%S13, &
805  glb%NumNP,elemstart,elemend,glb%NumElVol,glb%NumMatVol, &
806  glb%rho, glb%cd_fastest, glb%DetF_old, glb%VeloHalf, glb%Dt)
807 
808  ELSE
809 
810  CALL v3d10_nl(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,rnet,glb%Disp,glb%ci, &
811  glb%S11,glb%S22,glb%S33,glb%S12,glb%S23,glb%S13, &
812  glb%NumNP,elemstart,elemend,glb%NumElVol,glb%NumMatVol)
813  ENDIF
814  ELSE
815 
816  CALL v3d10r_nl(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,rnet,glb%Disp,glb%ci,glb%cj, &
817  glb%S11,glb%S22,glb%S33,glb%S12,glb%S23,glb%S13, &
818  glb%NumNP,elemstart,elemend,glb%NumElVol,glb%NumMatVol)
819  ENDIF
820 
821  ELSE IF (ianalysistype.EQ.2)THEN
822 
823 
824  IF(glb%iElIntgratn.EQ.0)THEN
825 !------- Linear elastic
826  CALL v3d10(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,rnet,glb%Disp,glb%ci, &
827  glb%S11,glb%S22,glb%S33,glb%S12,glb%S23,glb%S13, &
828  glb%NumNP,elemstart,elemend,glb%NumElVol,glb%NumMatVol)
829 
830  ELSE
831 
832  CALL v3d10_b_bar(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,rnet,glb%Disp,glb%ci, &
833  glb%S11,glb%S22,glb%S33,glb%S12,glb%S23,glb%S13, &
834  glb%NumNP,elemstart,elemend,glb%NumElVol,glb%NumMatVol)
835 
836  ENDIF
837  ENDIF
838  endif
839 
840  END SUBROUTINE internalforce_v3d10
841 
842  SUBROUTINE internalforce_v3d8( glb, Rnet, ElemStart, ElemEnd, iAnalysisType)
843 
844  IMPLICIT NONE
845 
846  TYPE(rocfrac_global) :: glb
847  REAL*8, DIMENSION(1:3*glb%NumNP) :: rnet
848 
849  INTEGER, INTENT(IN) :: elemstart, elemend, ianalysistype
850 
851 !------- Linear elastic
852  CALL v3d8_me(glb%MeshCoor,glb%MatIdVol,glb%ElConnVol,rnet,glb%Disp,glb%dmat, &
853  glb%S11,glb%S22,glb%S33,glb%S12,glb%S23,glb%S13, &
854  glb%NumNP,elemstart,elemend,glb%NumElVol,glb%NumMatVol,glb%Aenh,glb%enhanced_map,glb%mixed_map)
855 
856  END SUBROUTINE internalforce_v3d8
857 
858 END MODULE updatestructuralsoln
859 
subroutine v3d10(coor, matcstet, lmcstet, R_in, d, ci, S11l, S22l, S33l, S12l, S23l, S13l, numnp, nstart, nend, numcstet, numat_vol)
Definition: v3d10.f90:53
subroutine v3d4_damping(vhalf, R_in, numnp, numcstet, lmcstet, coor, nstart, nend, KappaDamp)
subroutine v3d10_nl(coor, matcstet, lmcstet, R_in, d, ci, S11, S22, S33, S12, S23, S13, numnp, nstart, nend, numcstet, numat_vol)
Definition: v3d10_nl.f90:53
subroutine c3d6nm(nsubn, nsubf, nsubn2, nsubf2, CoorOverlay, ElConnOverlay, sd_subface_parents, sd_subface_parents2, sd_subface_mate, sd_subface_nat_coors, sd_subface_nat_coors_mate, FaceOfVolEl, FaceOfVolEl2, NumNp, NumEl, ElConn, nf, nf2, MapFaceEl2Vol, MapFaceEl2Vol2, deltan, deltat, sigmax, taumax, Sinit, Rnet, Disp, Sthresh, NumMatCoh, MatID, SignFlag)
Definition: c3d6nm.f90:54
subroutine v3d4_nl_arruda_boyce(coor, matcstet, lmcstet, R_in, d, S11, S22, S33, S12, S23, S13, numnp, nstart, nend, numcstet, numat_vol, mu, kapp)
subroutine v3d10_b_bar(coor, matcstet, lmcstet, R_in, d, ci, S11l, S22l, S33l, S12l, S23l, S13l, numnp, nstart, nend, numcstet, numat_vol)
Definition: v3d10_B_Bar.f90:53
subroutine v3d10_nl_arruda_boyce(coor, matcstet, lmcstet, R_in, d, S11, S22, S33, S12, S23, S13, numnp, nstart, nend, numcstet, numat_vol, mu, kappa)
subroutine v3d10_thermal(NumEl, NumNP, ElConnVol, Coor, Kappa, Rnet, T, Rho, Cp, matcstet, numat_vol, MeshVel, ElemStart, ElemEnd)
j indices k indices k
Definition: Indexing.h:6
subroutine v3d10_nl_huang(coor, matcstet, lmcstet, R_in, d, S11, S22, S33, S12, S23, S13, numnp, nstart, nend, numcstet, numat_vol, STATEV_Part1, STATEV_Part2, NSTATEV, MATRIX, NMATRIX, PARTICLE, NPARTICLE, NPARTICLETYPE, INTERFAC, NINTERFAC, StrainTrace)
subroutine internalforce_v3d4(glb, Rnet, ElemStart, ElemEnd, iAnalysisType)
subroutine updatestructural(glb, NumProcs, Rnet)
subroutine v3d10_thermalexp(coor, matcstet, lmcstet, R_in, ci, S11l, S22l, S33l, numnp, nstart, nend, numcstet, numat_vol, CoeffExp, Temperature, Temperature0)
subroutine v3d4_nl(coor, matcstet, lmcstet, R_in, d, ci, S11, S22, S33, S12, S23, S13, numnp, nstart, nend, numcstet, numat_vol, xmu, xlambda)
Definition: v3d4_nl.f90:53
subroutine internalforce_v3d8(glb, Rnet, ElemStart, ElemEnd, iAnalysisType)
subroutine v3d10r_nl(coor, matcstet, lmcstet, R_in, d, ci, cj, S11, S22, S33, S12, S23, S13, numnp, nstart, nend, numcstet, numat_vol)
Definition: v3d10r_nl.f90:53
subroutine v3d4_thermalexp2(coor, matcstet, lmcstet, R_in, d, ci, S11, S22, S33, S12, S23, S13, numnp, nstart, nend, numcstet, numat_vol, CoeffExp, Temperature, Temperature0)
subroutine updatestructuralht(glb, NumProcs, Rnet, RnetHT)
subroutine v3d4_thermal(NumEl, NumNP, ElConnVol, Coor, Kappa, Rnet, T, Rho, Cp, matcstet, numat_vol, MeshVelo, ElemStart, ElemEnd)
subroutine v3d4_neohookeanincompress(coor, matcstet, lmcstet, R_in, d, S11, S22, S33, S12, S23, S13, numnp, nstart, nend, numcstet, numat_vol, xmu, xkappa)
blockLoc i
Definition: read.cpp:79
subroutine v3d10_nl_damping(coor, matcstet, lmcstet, R_in, d, ci, S11, S22, S33, S12, S23, S13, numnp, nstart, nend, numcstet, numat_vol, rho, cd_fastest, DetFold, velo, dt)
subroutine v3d4(coor, matcstet, lmcstet, R_in, d, ci, S11, S22, S33, S12, S23, S13, numnp, nstart, nend, numcstet, numat_vol)
Definition: v3d4.f90:53
subroutine v3d10_nl_matous(coor, lmcstet, R_in, d, S11, S22, S33, S12, S23, S13, numnp, nstart, nend, numcstet, p1, p2, Yin, SoftParam, Km, K2, Gm, G2, vm, nu2, a_eta, a_zeta, cm, c2, cd, Lo, Lm, L2, Mm, M2, L_bar, StrainOld)
virtual std::ostream & print(std::ostream &os) const
j indices j
Definition: Indexing.h:6
subroutine v3d8_me(coor, MatType, ElConnVol, Rnet, disp, dmat, S11, S22, S33, S12, S23, S13, numnp, nstart, nend, NumEL, NumMatType, Aenh, enhanced_map, mixed_map)
Definition: v3d8_me.f90:53
subroutine internalforce_v3d4ht(glb, Rnet, ElemStart, ElemEnd, iAnalysisType, Temperature)
subroutine v3d10_ale(v_bar, a_bar, d, vhalf, Rnet, E, xnu, rho, numnp, numat_vol, numlstet, matlstet, lmlstet, meshcoor, nstart, nend)
Definition: v3d10_ale.f90:65
subroutine internalforce_v3d10(glb, Rnet, ElemStart, ElemEnd, iAnalysisType)
subroutine v3d10_nl_arruda_boyce_damping(coor, matcstet, lmcstet, R_in, d, S11, S22, S33, S12, S23, S13, numnp, nstart, nend, numcstet, numat_vol, mu, kappa, rho, cd_fastest, DetFold, velo, dt)
subroutine v3d4_ale(v_bar, a_bar, d, vhalf, Rnet, E, xnu, rho, numnp, numat_vol, numcstet, matcstet, lmcstet, meshcoor, nstart, nend)
Definition: v3d4_ale.f90:53