Rocstar  1.0
Rocstar multiphysics simulation application
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
feminp.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 SUBROUTINE feminp(glb, myid)
54 
55  USE rocstar_rocfrac
56 
57 !!****f* Rocfrac/Source/feminp.f90
58 !!
59 !! NAME
60 !! feminp
61 !!
62 !! FUNCTION
63 !!
64 !! READ INPUT INFORMATION (i.e. Analysis Deck File)
65 !!
66 !! INPUTS
67 !! glb -- global array
68 !! myid -- processor id (starting at 0)
69 !!
70 !!****
71 
72  IMPLICIT NONE
73  include 'mpif.h'
74  TYPE(rocfrac_global) :: glb
75 
76 ! -- local
77  CHARACTER*200 :: keywd ! keyword parameter
78 
79  INTEGER :: ios ! io error
80  INTEGER :: ierr ! mpi error
81 
82  INTEGER :: i ! loop counter
83 ! -- global
84  INTEGER :: myid ! processor id
85  INTEGER :: lll ! length of keyword, with no trailing blanks
86 
87  REAL*8, ALLOCATABLE, DIMENSION(:) :: tmp_e, tmp_xnu, tmp_rho, tmp_alpha
88  INTEGER, ALLOCATABLE,DIMENSION(:) :: tmp_isolntype
89  REAL*8, ALLOCATABLE, DIMENSION(:) :: tmp_cp, tmp_kappaht
90  INTEGER, ALLOCATABLE,DIMENSION(:) :: tmp_isolntypeht
91 
92  ALLOCATE(tmp_e(1:10), tmp_xnu(1:10), tmp_rho(1:10), tmp_alpha(1:10), tmp_isolntype(1:10) )
93  ALLOCATE(tmp_cp(1:10), tmp_kappaht(1:10),tmp_isolntypeht(1:10))
94 !
95 ! -- Set Default values
96 
97 ! IF(myid.EQ.0) PRINT*,'ROCFRAC:: INSIDE feminp.f90'
98 
99  glb%DummyTractVal = 0.d0
100  glb%DummyBurnRate = 0.d0
101 
102  glb%NumNodeIO = 0
103  glb%ALEenabled = .false. ! ALE disabled by default
104  glb%ipstatic = .false. ! don't subract the intial pressure
105  glb%ReStart = .false.
106  glb%iElType = 4 ! default to 4 node tets
107  glb%IONEWER = .false.
108  glb%DampEnabled = .false.
109  glb%iElIntgratn = 0
110  glb%NdMassLump = 0
111  glb%EnforceTractionS = .false.! To enforce traction where no fluids pressure
112  glb%NumEntries = 0
113  glb%NumMatVol = 0
114  glb%NumMatVolHT = 0
115  glb%NumMatCoh = 0
116  glb%cd_fastest = 0.d0
117  glb%NdBasedEl = .false.
118  glb%UnDefConfig = .false.
119  glb%HeatTransSoln = .false.
120  glb%Temperature0 = 0.
121  glb%ArtificialDamping = .false.
122  glb%EnforceTractionS = .false.
123  glb%EnforceTractionSF = .false.
124  glb%DebondPart = .false.
125  glb%DebondPart_MATOUS = .false.
126  glb%ThermalExpansion = .false.
127  glb%AmplitudeTable = .false.
128  glb%debug_state = .false.
129  glb%NumProbesEl = 0
130  glb%NumProbesNd = 0
131  glb%OverlayExist = .false.
132  glb%DummyFlux = 0.d0
133 
134  glb%NumElOverlay = 0
135  glb%NumNpOverlay = 0
136  glb%Verb = 1
137 !
138 ! -- Open Analysis Deck File
139 
140  OPEN(glb%io_input,file='./Rocfrac/RocfracControl.txt',status='old',iostat=ios)
141  IF(ios.NE.0)THEN
142  IF(myid.EQ.0) print*, 'ROCFRAC:: Unable to find RocfracControl.txt - STOPPING'
143  CALL mpi_finalize(glb%MPI_COMM_ROCFRAC,ierr)
144  stop
145  ENDIF
146 !
147 ! Input deck summary file
148 
149  IF(myid.EQ.0)THEN
150  OPEN(glb%io_sum,file='Rocfrac/Modin/InputSummary.res',status='unknown',iostat=ios)
151  IF(ios.NE.0)THEN ! Try without Modin
152  OPEN(glb%io_sum,file='Rocfrac/InputSummary.res',status='unknown',iostat=ios)
153  END IF
154 
155  IF(ios.NE.0)THEN
156  IF(myid.EQ.0) print*, 'ROCFRAC:: Unable to find InputSummary.res under Rocfrac/Modin/ or Rocfrac/ - STOPPING'
157  CALL mpi_finalize(glb%MPI_COMM_ROCFRAC,ierr)
158  stop
159  ENDIF
160  ENDIF
161 
162 !
163 ! -- Read Analysis Deck File
164 
165  rewind glb%io_input
166 1 READ(glb%io_input,'(A)',iostat=ios) keywd
167 ! if(myid.EQ.0) print*,keywd
168  IF(ios.LT.0) THEN ! Negative ios means end-of-file
169  print*,' *END parameter not found - STOPPING'
170  CALL mpi_finalize(glb%MPI_COMM_ROCFRAC,ierr)
171  stop
172  ENDIF
173 
174 2 CONTINUE
175 
176 !
177 ! Comment field
178 
179  IF(keywd(1:1).NE.'*')THEN
180  goto 1
181  ENDIF
182 
183  IF(keywd(1:2).EQ.'**')THEN
184  goto 1
185  ENDIF
186 
187  lll = len_trim(keywd)
188 
189  IF(myid.EQ.0) WRITE(glb%io_sum,'(2X,A,A)') keywd(1:lll)
190 
191  IF(keywd(1:4).EQ.'*END') THEN
192  goto 3
193 
194 ! Amplitude
195 
196  ELSE IF(keywd(1:10).EQ.'*AMPLITUDE') THEN
197  CALL amplitude_sub(glb,keywd)
198  goto 1
199 
200  ELSE IF(keywd(1:13).EQ.'*HYPERELASTIC')THEN
201  CALL matmodel_hyperelastic(glb, keywd, &
202  tmp_e, tmp_xnu, tmp_rho, tmp_alpha, tmp_isolntype)
203  goto 1
204 
205  ELSE IF(keywd(1:14).EQ.'*HEAT TRANSFER')THEN
206  CALL heat_transfer_sub(glb, keywd, tmp_kappaht, tmp_cp, tmp_isolntypeht)
207  goto 1
208 
209  ELSE IF(keywd(1:8).EQ.'*ELASTIC')THEN
210  CALL matmodel_elastic(glb, keywd, &
211  tmp_e, tmp_xnu, tmp_rho, tmp_alpha, tmp_isolntype)
212  goto 1
213 
214  ELSE IF(keywd(1:8).EQ.'*ELEMENT')THEN
215  CALL element_sub(glb,keywd)
216  goto 1
217  ELSE IF(keywd(1:8).EQ.'*DAMPING') THEN
218  READ(glb%io_input,*) glb%KappaDamp
219  glb%DampEnabled = .true.
220  goto 1
221 !
222 ! Prefix for input/output files
223 !
224  ELSE IF(keywd(1:7).EQ.'*PREFIX') THEN
225  CALL prefix_sub(glb)
226  goto 1
227  ELSE IF(keywd(1:4).EQ.'*ALE') THEN
228  CALL ale_sub(glb)
229  goto 1
230  ELSE IF(keywd(1:5).EQ.'*NRUN') THEN
231  CALL nrun_sub(glb)
232  goto 1
233  ELSE IF(keywd(1:8).EQ.'*DYNAMIC') THEN
234  CALL dynamic_sub(glb,keywd)
235  goto 1
236 !*** OBSOLETE REMOVE
237  ELSE IF(keywd(1:7).EQ.'*MATVOL') THEN ! Read volumetric material props.
238  CALL matvol_sub(glb, &
239  tmp_e, tmp_xnu, tmp_rho, tmp_alpha, tmp_isolntype)
240  goto 1
241  ELSE IF(keywd(1:7).EQ.'*MATCOH') THEN ! Read cohesive material props
242  CALL matcoh_sub(glb)
243  goto 1
244  ELSE IF(keywd(1:7).EQ.'*PLOAD1') THEN
245  CALL pload1_sub(glb)
246  goto 1
247  ELSE IF(keywd(1:5).EQ.'*NODE') THEN
248  CALL nodeio_sub(myid)
249  goto 1
250  ELSE IF(keywd(1:11).EQ.'*DUMMYTRACT') THEN
251  CALL dummytract_sub(glb,keywd)
252  goto 1
253  ELSE IF(keywd(1:10).EQ.'*DUMMYBURN') THEN
254  READ(glb%io_input,*) glb%DummyBurnRate
255  goto 1
256  ELSE IF(keywd(1:10).EQ.'*DUMMYFLUX') THEN
257  READ(glb%io_input,*) glb%DummyFlux
258  goto 1
259  ELSE IF(keywd(1:11).EQ.'*BOUNDARYMM') THEN
260  CALL boundarymm_sub(glb)
261  goto 1
262  ELSE IF(keywd(1:11).EQ.'*BOUNDARYHT') THEN
263  CALL boundaryht_sub(glb)
264  goto 1
265  ELSE IF(keywd(1:9).EQ.'*BOUNDARY') THEN
266  CALL boundary_sub(glb)
267  goto 1
268  ELSE IF(keywd(1:9).EQ.'*IPSTATIC')THEN ! to subract out the initial pressure
269  glb%ipstatic = .true.
270  goto 1
271  ELSE IF(keywd(1:8).EQ.'*IONEWER')THEN ! new gen 2.5 input format
272  glb%IONEWER = .true.
273  goto 1
274  ELSE IF(keywd(1:10).EQ.'*DEFCONFIG')THEN
275  glb%UnDefConfig = .true.
276  goto 1
277 
278  ELSE IF(keywd(1:16).EQ.'*MICROMECHANICAL') THEN ! micromechanical model
279  CALL micromechanical_sub(glb, keywd, &
280  tmp_e, tmp_xnu, tmp_rho, tmp_alpha, tmp_isolntype)
281  goto 1
282 
283  ELSE IF(keywd(1:11).EQ.'*ARTDAMPING')THEN
284  glb%ArtificialDamping = .true.
285  goto 1
286  ELSE IF(keywd(1:6).EQ.'*PROBE')THEN
287  CALL probe_sub(glb)
288  goto 1
289  ELSE IF(keywd(1:18).EQ.'*INITIAL CONDITION')THEN
290  CALL initialcondition_sub(glb,keywd)
291  goto 1
292  ELSE IF(keywd(1: ).EQ.'*VERBOSITY')THEN
293  READ(glb%io_input,*) glb%Verb
294  goto 1
295  ELSE IF(keywd(1:6).EQ.'*DEBUG')THEN
296  glb%debug_state = .true.
297  goto 1
298 
299  ELSE IF(keywd(1:1).EQ.'*')THEN
300  print*,'ROCFRAC: ERROR'
301  print*,'ROCFRAC: CONTROL DECK OPTION ',trim(keywd),' NOT SUPPORTED'
302  print*,'ROCFRAC: STOPPING'
303  CALL mpi_barrier(glb%MPI_COMM_ROCFRAC,ierr)
304  CALL mpi_finalize(glb%MPI_COMM_ROCFRAC,ierr)
305  stop
306  ELSE
307  goto 1
308  ENDIF
309 3 CONTINUE
310 
311  CLOSE(glb%io_input)
312 !
313 !---- Write analysis summary file
314 
315  ALLOCATE(glb%E (1:glb%NumMatVol) )
316  ALLOCATE(glb%xnu (1:glb%NumMatVol) )
317  ALLOCATE(glb%rho (1:glb%NumMatVol) )
318  ALLOCATE(glb%alpha (1:glb%NumMatVol) )
319  ALLOCATE(glb%xmu (1:glb%NumMatVol) )
320  ALLOCATE(glb%xlambda (1:glb%NumMatVol) )
321  ALLOCATE(glb%xkappa (1:glb%NumMatVol) )
322 
323  ALLOCATE(glb%iSolnType(1:glb%NumMatVol) )
324 
325 
326  IF(glb%DebondPart_Matous) glb%rho(1) = tmp_rho(1)
327 
328  DO i = 1, glb%NumMatVol
329  glb%E(i) = tmp_e(i)
330  glb%xnu(i) = tmp_xnu(i)
331  glb%rho(i) = tmp_rho(i)
332  glb%alpha(i) = tmp_alpha(i)
333  glb%iSolnType(i) = tmp_isolntype(i)
334 
335  glb%xmu(i) = glb%E(i)/(2.d0*(1.d0+glb%xnu(i)))
336 
337  glb%xlambda(i) = 2.d0*glb%xmu(i)*glb%xnu(i)/(1.d0-2.d0*glb%xnu(i))
338 
339 ! bulk modulus
340  glb%xkappa(i) = glb%xlambda(i) + 2.d0/3.d0*glb%xmu(i)
341 
342  ENDDO
343 
344  IF(glb%HeatTransSoln)THEN
345  glb%ThermalDiffusivity = -1.
346 
347  ALLOCATE(glb%iSolnTypeHT(1:glb%NumMatVolHT) )
348  ALLOCATE(glb%KappaHT (1:glb%NumMatVolHT) )
349  ALLOCATE(glb%Cp (1:glb%NumMatVolHT) )
350 
351  DO i = 1, glb%NumMatVolHT
352  glb%KappaHT(i) = tmp_kappaht(i)
353  glb%Cp(i) = tmp_cp(i)
354  glb%iSolnTypeHT(i) = tmp_isolntypeht(i)
355 
356 
357 
358  glb%ThermalDiffusivity = max( glb%ThermalDiffusivity, glb%KappaHT(i)/(glb%rho(i)*glb%Cp(i)) )
359 
360 
361  ENDDO
362  ENDIF
363 
364  DEALLOCATE(tmp_kappaht, tmp_cp, tmp_isolntypeht)
365  DEALLOCATE(tmp_e,tmp_xnu,tmp_rho,tmp_alpha,tmp_isolntype)
366 
367  IF(myid.EQ.0)THEN
368  WRITE(glb%io_sum,50) glb%prefx
369  WRITE(glb%io_sum,60)
370  IF (glb%restart) WRITE(glb%io_sum,65)
371  WRITE(glb%io_sum,80) glb%NumMatVol,glb%NumMatCoh
372  WRITE(glb%io_sum,85) glb%CourantRatio
373  WRITE(glb%io_sum,86)
374  DO i = 1,glb%NumMatCoh
375  WRITE(glb%io_sum,88) glb%deltan(i),glb%deltat(i),glb%SigmaMax(i),glb%TauMax(i),glb%Sinit(i)
376 
377 
378  ENDDO
379  WRITE(glb%io_sum,90)
380  WRITE(glb%io_sum,91)
381  DO i = 1,glb%NumMatVol
382  WRITE(glb%io_sum,92) i,glb%E(i),glb%xnu(i),glb%rho(i),glb%alpha(i)
383  IF(glb%iSolnType(i).EQ.0)THEN
384  WRITE(glb%io_sum,*) ' Material Model = Arruda-Boyce '
385  ELSE IF(glb%iSolnType(i).EQ.-1)THEN
386  WRITE(glb%io_sum,*) ' Material Model = NeoHookean Incompressible '
387  ELSE IF(glb%iSolnType(i).EQ.1)THEN
388  WRITE(glb%io_sum,*) ' Material Model = Elastic, Large Deformation '
389  ELSE IF(glb%iSolnType(i).EQ.2)THEN
390  WRITE(glb%io_sum,*) ' Material Model = Elastic, Small Deformation '
391  ELSE
392  WRITE(glb%io_sum,*) ' Error, Not a valid Material model, = ', glb%iSolnType(i)
393  ENDIF
394 
395  ENDDO
396 
397  IF(glb%HeatTransSoln) WRITE(glb%io_sum,*) 'Heat Transfer Solution'
398  CLOSE(glb%io_sum)
399  ENDIF
400 
401  RETURN
402 
403 !--------------------------------FORMATS--------------------------------
404 !
405 
406 50 FORMAT(//,'DYNAMIC 3D LINEAR ELASTIC FEA',///,'Job id: ',a20,//)
407 60 FORMAT('*** ISOTROPIC ANALYSIS ***',/)
408 65 FORMAT(' *** RESTART ***',/)
409 80 FORMAT(1x,'Number of material types (NUMAT_VOL) =',i12 &
410  /,1x,'Number of material types (NUMAT_COH) =',i12)
411 85 FORMAT(/,1x,'Steps per characteristic length (STEPS) =',e12.4)
412 86 FORMAT(///,1x,'COHESIVE ELEMENT DATA')
413 88 FORMAT(/,1x,'Characteristic lengths for the cohesive law:', &
414  /,1x,' normal (DELTAN) =',e12.4, &
415  /,1x,' tangential (DELTAT) =',e12.4, &
416  /,1x,'Maximum normal stress (GLB%SIGMAMAX) =',e12.4, &
417  /,1x,'Maximum shearing stress (TAUMAX) =',e12.4, &
418  /,1x,'The initial Sthreshold S(init) =',e12.4)
419 90 FORMAT(///,1x,'VOLUMETRIC ELEMENT DATA')
420 91 FORMAT(/,4x,'MATERIAL SETS')
421 92 FORMAT(/,9x,'SET',4x,i4 &
422  /,12x,'E =',10x,e13.5, &
423  /,12x,'Nu =',10x,e13.5, &
424  /,12x,'rho =',10x,e13.5, &
425  /,12x,'alpha =',10x,e13.5)
426 
427 END SUBROUTINE feminp
428 
429 
430 
431 SUBROUTINE prefix_sub(glb)
432 
433 !!****f* Rocfrac/Source/feminp/PREFIX_SUB
434 !!
435 !! NAME
436 !! PREFIX_SUB
437 !!
438 !! FUNCTION
439 !!
440 !! Reads prefix keyword (i.e. Analysis Deck File)
441 !!
442 !! INPUTS
443 !! glb -- global array
444 !!
445 !!****
446 
447  USE rocstar_rocfrac
448 
449  IMPLICIT NONE
450 
451  TYPE(rocfrac_global) :: glb
452 
453  READ(glb%io_input,'(A)') glb%prefx
454  glb%prefx_lngth = len_trim(glb%prefx)
455 
456  RETURN
457 END SUBROUTINE prefix_sub
458 
459 SUBROUTINE ale_sub(glb)
460 
461 !!****f* Rocfrac/Source/feminp/ALE_SUB
462 !!
463 !! NAME
464 !! ALE_SUB
465 !!
466 !! FUNCTION
467 !!
468 !! ALE keyword turns on ALE routines
469 !!
470 !! INPUTS
471 !! glb -- global array
472 !!
473 !!****
474 
475  USE rocstar_rocfrac
476 
477  IMPLICIT NONE
478 
479  TYPE(rocfrac_global) :: glb
480 
481 ! kappa parameter for mesh motion
482 
483  READ(glb%io_input,*) glb%kappa
484  glb%ALEenabled = .true.
485 
486  !print*,'Running with ALE'
487 
488  RETURN
489 END SUBROUTINE ale_sub
490 
491 SUBROUTINE nrun_sub(glb)
492 
493  USE rocstar_rocfrac
494 
495  IMPLICIT NONE
496 
497  TYPE(rocfrac_global) :: glb
498 
499  READ(glb%io_input,*) glb%CourantRatio
500 
501 !RAF Make this backward compatible with old input decks --
502 !RAF You would not want a ratio > 1, right?
503 
504  IF (glb%CourantRatio > 1.0d0) THEN
505  glb%CourantRatio = 1.0d0 / glb%CourantRatio
506  ENDIF
507 
508  RETURN
509 END SUBROUTINE nrun_sub
510 
511 SUBROUTINE dynamic_sub(glb,keywd)
512 
513 !!****f* Rocfrac/Source/feminp/DYNAMIC_SUB
514 !!
515 !! NAME
516 !! DYNAMIC_SUB
517 !!
518 !! FUNCTION
519 !!
520 !! Courant limit multiplier
521 !!
522 !! INPUTS
523 !! glb -- global array
524 !! keywd -- keywd for control deck
525 !!
526 !!****
527 
528  USE rocstar_rocfrac
529 
530  IMPLICIT NONE
531 
532  CHARACTER(len=200) :: keywd
533  CHARACTER(len=16) :: scalefactor
534  INTEGER :: k1, k2
535 
536 
537  TYPE(rocfrac_global) :: glb
538 
539  CALL locchr(keywd,'SCALE FACTOR ',12,8,k1,k2)
540 
541  scalefactor = keywd(k1:k2)
542 
543  CALL rchar(scalefactor,glb%CourantRatio)
544 
545  RETURN
546 END SUBROUTINE dynamic_sub
547 
548 SUBROUTINE dummytract_sub(glb,keywd)
549 
550 !!****f* Rocfrac/Source/feminp/DUMMYTRACT_SUB
551 !!
552 !! NAME
553 !! DUMMYTRACT_SUB
554 !!
555 !! FUNCTION
556 !! Marks if applying the dummy traction to what surface
557 !!
558 !!
559 !! INPUTS
560 !! glb -- global array
561 !! keywd -- keywd for control deck
562 !!
563 !!****
564 
565  USE rocstar_rocfrac
566 
567  IMPLICIT NONE
568 
569  CHARACTER(len=200) :: keywd
570  INTEGER :: k1, k2
571  CHARACTER :: option*16
572 
573  TYPE(rocfrac_global) :: glb
574 
575  CALL locchr(keywd,'INTERFACE ',9,8,k1,k2)
576 
577  option = keywd(k1:k2)
578 
579  IF(option.EQ.'S') glb%EnforceTractionS = .true.
580  IF(option.EQ.'SF') glb%EnforceTractionSF = .true.
581  IF(option.EQ.'ALL')THEN
582  glb%EnforceTractionS = .true.
583  glb%EnforceTractionSF = .true.
584  ENDIF
585 
586  READ(glb%io_input,*) glb%DummyTractVal
587 
588  RETURN
589 END SUBROUTINE dummytract_sub
590 
591 SUBROUTINE matmodel_hyperelastic(glb,keywd, &
592  tmp_e, tmp_xnu, tmp_rho, tmp_alpha, tmp_isolntype)
593 
594 !!****f* Rocfrac/Source/feminp/MATMODEL_HYPERELASTIC
595 !!
596 !! NAME
597 !! MATMODEL_HYPERELASTIC
598 !!
599 !! FUNCTION
600 !! reads material model type, and the
601 !! hyperelastic material parameters
602 !!
603 !!
604 !! INPUTS
605 !! glb -- global array
606 !! keywd -- keywd for control deck
607 !!
608 !! OUTPUT
609 !! tmp_E -- Young's Modulus
610 !! tmp_xnu -- Possion's ratio
611 !! tmp_rho -- Density
612 !! tmp_alpha -- thermal coefficient of expansion
613 !! tmp_iSolnType -- material type model
614 !!
615 !!****
616 
617  USE rocstar_rocfrac
618 
619  IMPLICIT NONE
620 
621  TYPE(rocfrac_global) :: glb
622 
623  CHARACTER :: keywd*200
624  INTEGER :: key ! 0 = no key word found, 1 = key word found)
625  CHARACTER*26 :: mattype
626  INTEGER :: nummattype
627  INTEGER :: i,ii,ierr
628  REAL*8, DIMENSION(1:10) :: tmp_e, tmp_xnu, tmp_rho, tmp_alpha
629  INTEGER, DIMENSION(1:10) :: tmp_isolntype
630 
631  CALL conchr(keywd,'NEOHOOKINC ',10,13,key)
632  IF(key.EQ.1) mattype = 'NEOHOOKINC'
633 
634  CALL conchr(keywd,'ARRUDA-BOYCE ',12,13,key)
635  IF(key.EQ.1) mattype = 'ARRUDA-BOYCE'
636 
637  SELECT CASE (trim(mattype))
638 
639  CASE ('ARRUDA-BOYCE')
640 
641  READ(glb%io_input,*) nummattype
642 
643  DO i = 1, nummattype
644  glb%NumMatVol = glb%NumMatVol + 1
645 
646  ii = glb%NumMatVol
647 
648  IF(ii.GT.10)THEN
649  print*,'ROCFRAC :: ERROR'
650  print*,'Number of materials GREATER then 10'
651  CALL mpi_finalize(glb%MPI_COMM_ROCFRAC,ierr)
652  ENDIF
653 
654  READ(glb%io_input,*) tmp_e(ii), tmp_xnu(ii), tmp_rho(ii), tmp_alpha(ii)
655 
656  glb%cd_fastest = max( glb%cd_fastest, &
657  sqrt(tmp_e(ii)*(1.d0-tmp_xnu(ii))/tmp_rho(ii)/(1.d0+tmp_xnu(ii))/(1.d0-2.d0*tmp_xnu(ii) )) )
658  tmp_isolntype(ii) = 0
659 
660  ENDDO
661 
662  CASE ('NEOHOOKINC')
663 
664  READ(glb%io_input,*) nummattype
665 
666  DO i = 1, nummattype
667  glb%NumMatVol = glb%NumMatVol + 1
668 
669  ii = glb%NumMatVol
670 
671  IF(ii.GT.10)THEN
672  print*,'ROCFRAC :: ERROR'
673  print*,'Number of materials GREATER then 10'
674  CALL mpi_finalize(glb%MPI_COMM_ROCFRAC,ierr)
675  ENDIF
676 
677  READ(glb%io_input,*) tmp_e(ii), tmp_xnu(ii), tmp_rho(ii), tmp_alpha(ii)
678 
679  glb%cd_fastest = max( glb%cd_fastest, &
680  sqrt(tmp_e(ii)*(1.d0-tmp_xnu(ii))/tmp_rho(ii)/(1.d0+tmp_xnu(ii))/(1.d0-2.d0*tmp_xnu(ii) )) )
681  tmp_isolntype(ii) = -1
682 
683  ENDDO
684 
685  CASE default
686 
687  print*,'ROCFRAC :: ERROR'
688  print*,'ROCFRAC :: *HYPERELASTIC KEYWORD ',trim(mattype), ' NOT FOUND'
689  print*,'ROCFRAC :: STOPPING'
690  CALL mpi_finalize(glb%MPI_COMM_ROCFRAC,ierr)
691 
692 
693  END SELECT
694 
695 
696 END SUBROUTINE matmodel_hyperelastic
697 
698 SUBROUTINE matmodel_elastic(glb,keywd, &
699  tmp_e, tmp_xnu, tmp_rho, tmp_alpha, tmp_isolntype)
700 
701 !!****f* Rocfrac/Source/feminp/MATMODEL_ELASTIC
702 !!
703 !! NAME
704 !! MATMODEL_ELASTIC
705 !!
706 !! FUNCTION
707 !! reads material model type, and the
708 !! elastic material parameters
709 !!
710 !!
711 !! INPUTS
712 !! glb -- global array
713 !! keywd -- keywd for control deck
714 !!
715 !! OUTPUT
716 !! tmp_E -- Young's Modulus
717 !! tmp_xnu -- Possion's ratio
718 !! tmp_rho -- Density
719 !! tmp_alpha -- thermal coefficient of expansion
720 !! tmp_iSolnType -- material type model
721 !!
722 !!****
723 
724  USE rocstar_rocfrac
725 
726  IMPLICIT NONE
727 
728  TYPE(rocfrac_global) :: glb
729 
730  CHARACTER :: keywd*200
731  INTEGER :: key ! 0 = no key word found, 1 = key word found)
732  CHARACTER*26 :: mattype
733  INTEGER :: nummattype
734  INTEGER :: i,ii,ierr
735  REAL*8, DIMENSION(1:10) :: tmp_e, tmp_xnu, tmp_rho, tmp_alpha
736  INTEGER, DIMENSION(1:10) :: tmp_isolntype
737 
738  INTEGER :: k1, k2
739  CHARACTER :: nlgeom*16
740  INTEGER :: nlgeomtype
741 
742  CALL locchr(keywd,'NLGEOM ',6,8,k1,k2)
743 
744  nlgeom = keywd(k1:k2)
745 
746  nlgeomtype = 1 ! defaults to nl
747  IF(nlgeom.EQ.'NO')THEN
748  nlgeomtype = 2
749  ENDIF
750 
751 
752  READ(glb%io_input,*) nummattype
753 
754  DO i = 1, nummattype
755 
756  glb%NumMatVol = glb%NumMatVol + 1
757 
758  ii = glb%NumMatVol
759 
760  IF(ii.GT.10)THEN
761  print*,'ROCFRAC :: ERROR'
762  print*,'Number of ELASTIC materials GREATER then 10'
763  CALL mpi_finalize(glb%MPI_COMM_ROCFRAC,ierr)
764  ENDIF
765 
766  READ(glb%io_input,*) tmp_e(ii), tmp_xnu(ii), tmp_rho(ii), tmp_alpha(ii)
767 
768  IF(tmp_alpha(ii).NE.0.d0) glb%ThermalExpansion = .true.
769 
770  glb%cd_fastest = max( glb%cd_fastest, &
771  sqrt(tmp_e(ii)*(1.d0-tmp_xnu(ii))/tmp_rho(ii)/(1.d0+tmp_xnu(ii))/(1.d0-2.d0*tmp_xnu(ii) )) )
772  tmp_isolntype(ii) = nlgeomtype
773 
774  ENDDO
775 
776 END SUBROUTINE matmodel_elastic
777 
778 !*** OBSOLETE REMOVE
779 
780 SUBROUTINE matvol_sub(glb, &
781  tmp_e, tmp_xnu, tmp_rho, tmp_alpha, tmp_isolntype)
782 
783  USE rocstar_rocfrac
784 
785  IMPLICIT NONE
786 
787  TYPE(rocfrac_global) :: glb
788 
789  INTEGER :: i,ii,ierr ! loop counter
790  integer :: nummattype
791  REAL*8, DIMENSION(1:10) :: tmp_e, tmp_xnu, tmp_rho, tmp_alpha
792  INTEGER, DIMENSION(1:10) :: tmp_isolntype
793 
794  READ(glb%io_input,*) nummattype
795 
796  DO i = 1, nummattype
797  glb%NumMatVol = glb%NumMatVol + 1
798 
799  ii = glb%NumMatVol
800 
801  IF(ii.GT.10)THEN
802  print*,'ROCFRAC :: ERROR'
803  print*,'Number of materials GREATER then 10'
804  CALL mpi_finalize(glb%MPI_COMM_ROCFRAC,ierr)
805  stop
806  ENDIF
807 
808  READ(glb%io_input,*) tmp_e(ii), tmp_xnu(ii), tmp_rho(ii), tmp_alpha(ii),tmp_isolntype(ii)
809 
810 
811  glb%cd_fastest = max( glb%cd_fastest, &
812  sqrt(tmp_e(ii)*(1.d0-tmp_xnu(ii))/tmp_rho(ii)/(1.d0+tmp_xnu(ii))/(1.d0-2.d0*tmp_xnu(ii) )) )
813  ENDDO
814  RETURN
815 END SUBROUTINE matvol_sub
816 
817 SUBROUTINE matcoh_sub(glb)
818 
819  USE rocstar_rocfrac
820 
821  IMPLICIT NONE
822 
823  TYPE(rocfrac_global) :: glb
824 
825  INTEGER :: i ! loop counter
826 
827  READ(glb%io_input,*) glb%NumMatCoh
828 
829  ALLOCATE(glb%deltan(1:glb%NumMatCoh))
830  ALLOCATE(glb%deltat(1:glb%NumMatCoh))
831  ALLOCATE(glb%SigmaMax(1:glb%NumMatCoh))
832  ALLOCATE(glb%TauMax(1:glb%NumMatCoh))
833  ALLOCATE(glb%Sinit (1:glb%NumMatCoh))
834 
835  DO i = 1, glb%NumMatCoh
836  READ(glb%io_input,*) glb%deltan(i),glb%deltat(i),glb%SigmaMax(i),glb%TauMax(i), glb%Sinit(i)
837  ENDDO
838 
839  RETURN
840 END SUBROUTINE matcoh_sub
841 
842 SUBROUTINE pload1_sub(glb)
843 
844  USE rocstar_rocfrac
845 
846  IMPLICIT NONE
847 
848  TYPE(rocfrac_global):: glb
849 
850  INTEGER :: pload1_input ! fix not donei ! loop counter
851 
852  READ(glb%io_input,*) pload1_input
853 
854  RETURN
855 END SUBROUTINE pload1_sub
856 
857 SUBROUTINE nodeio_sub(myid)
858 
859  USE rocstar_rocfrac
860 
861  IMPLICIT NONE
862 
863  TYPE(rocfrac_global) :: glb
864 
865  INTEGER :: i ! loop counter
866  INTEGER :: myid
867 
868  READ(glb%io_input,*) glb%NumNodeIO, glb%NumNodeIOpid
869 
870  IF(glb%NumNodeIOpid.EQ.myid)THEN
871 
872  ALLOCATE(glb%NodeIO(1:glb%NumNodeIO))
873 
874  DO i = 1, glb%NumNodeIO
875  READ(glb%io_input,*) glb%NodeIO
876  ENDDO
877 
878  ELSE
879 
880  glb%NumNodeIO = 0
881 
882  ENDIF
883 
884  RETURN
885 END SUBROUTINE nodeio_sub
886 
887 
888 SUBROUTINE boundary_sub(glb)
889 
890 !!****f* Rocfrac/Source/feminp/BOUNDARY_SUB
891 !!
892 !! NAME
893 !! BOUNDARY_SUB
894 !!
895 !! FUNCTION
896 !! This option is used to prescibe boundary conditions
897 !! at nodes.
898 !!
899 !!
900 !! INPUTS
901 !! glb -- global array
902 !!
903 !!****
904 
905  USE rocstar_rocfrac
906 
907  IMPLICIT NONE
908 
909  TYPE(rocfrac_global) :: glb
910 
911  INTEGER :: i,iaux ! loop counter
912  INTEGER :: numbcflags
913 
914  READ(glb%io_input,*) numbcflags
915 
916  ALLOCATE(glb%bcCond(1:numbcflags))
917 
918  DO i = 1, numbcflags
919  READ(glb%io_input,*) iaux,glb%bcCond(i)%BCtypeX,glb%bcCond(i)%BCtypeY,glb%bcCond(i)%BCtypeZ, &
920  glb%bcCond(i)%BCvalueX,glb%bcCond(i)%BCvalueY,glb%bcCond(i)%BCvalueZ
921  ENDDO
922 
923  RETURN
924 END SUBROUTINE boundary_sub
925 
926 SUBROUTINE boundaryht_sub(glb)
927 
928 !!****f* Rocfrac/Source/feminp/BOUNDARYHT_SUB
929 !!
930 !! NAME
931 !! BOUNDARYHT_SUB
932 !!
933 !! FUNCTION
934 !! This option is used to prescibe heat transfer
935 !! boundary conditions at nodes.
936 !!
937 !!
938 !! INPUTS
939 !! glb -- global array
940 !!
941 !!****
942 
943  USE rocstar_rocfrac
944 
945  IMPLICIT NONE
946 
947  TYPE(rocfrac_global) :: glb
948 
949  INTEGER :: i,iaux ! loop counter
950  INTEGER :: numbcflagsht
951 
952  READ(glb%io_input,*) numbcflagsht
953 
954 
955  ALLOCATE(glb%bcCondHT(1:numbcflagsht))
956 
957  DO i = 1, numbcflagsht
958  READ(glb%io_input,*) iaux,glb%bcCondHT(i)%BCtypeX,glb%bcCondHT(i)%BCtypeY,glb%bcCondHT(i)%BCtypeZ, &
959  glb%bcCondHT(i)%BCvalueX,glb%bcCondHT(i)%BCvalueY,glb%bcCondHT(i)%BCvalueZ
960  ENDDO
961 
962  RETURN
963 END SUBROUTINE boundaryht_sub
964 
965 SUBROUTINE boundarymm_sub(glb)
966 
967 !!****f* Rocfrac/Source/feminp/BOUNDARYMM_SUB
968 !!
969 !! NAME
970 !! BOUNDARYMM_SUB
971 !!
972 !! FUNCTION
973 !! This option is used to prescibe mesh motion boundary
974 !! conditions at nodes.
975 !!
976 !! INPUTS
977 !! glb -- global array
978 !!
979 !!****
980 
981  USE rocstar_rocfrac
982 
983  IMPLICIT NONE
984 
985  TYPE(rocfrac_global) :: glb
986 
987  INTEGER :: i,iaux ! loop counter
988  INTEGER :: numbcflagsmm
989 
990  READ(glb%io_input,*) numbcflagsmm
991 
992  ALLOCATE(glb%bcCondmm(1:numbcflagsmm))
993 
994  DO i = 1, numbcflagsmm
995  READ(glb%io_input,*) iaux,glb%bcCondmm(i)%BCtypeX,glb%bcCondmm(i)%BCtypeY,glb%bcCondmm(i)%BCtypeZ, &
996  glb%bcCondmm(i)%BCvalueX,glb%bcCondmm(i)%BCvalueY,glb%bcCondmm(i)%BCvalueZ
997  ENDDO
998 
999  RETURN
1000 END SUBROUTINE boundarymm_sub
1001 
1002 SUBROUTINE amplitude_sub(glb,keywd)
1003 
1004  USE rocstar_rocfrac
1005 
1006  IMPLICIT NONE
1007 
1008  TYPE(rocfrac_global) :: glb
1009 
1010  CHARACTER :: keywd*200
1011  INTEGER :: k1, k2
1012 
1013  INTEGER :: i ! loop counter
1014 
1015  CHARACTER :: amptype*16
1016  REAL*8 :: der1
1017  REAL*8 :: slp, intcpt, rise, run
1018  REAL*8, POINTER, DIMENSION(:,:) :: tableval
1019  REAL*8 :: x1, y1, x2, y2
1020 
1021  CALL locchr(keywd,'DEFINITION ',10,10,k1,k2)
1022 
1023  amptype = keywd(k1:k2)
1024 
1025  IF(amptype.EQ.'TABULAR')THEN
1026 
1027  glb%AmplitudeTable = .true.
1028 
1029  READ(glb%io_input,*) glb%NumEntries
1030 
1031 ! Time
1032 ! 1st Derivative
1033 ! 2nd Derivative
1034 
1035  ALLOCATE(tableval(1:2,glb%NumEntries))
1036  DO i = 1, glb%NumEntries
1037  READ(glb%io_input,*) tableval(1:2,i)
1038  ENDDO
1039 
1040  glb%NumEntries = glb%NumEntries - 1
1041 
1042  ALLOCATE( glb%AmpTable(1:3,glb%NumEntries))
1043 
1044  DO i = 1, glb%NumEntries
1045 
1046  x1 = tableval(1,i)
1047  y1 = tableval(2,i)
1048  x2 = tableval(1,i+1)
1049  y2 = tableval(2,i+1)
1050 
1051  rise = y2-y1
1052  run = x2-x1
1053 
1054  slp = (y2-y1)/(x2-x1)
1055 
1056  IF(abs(slp*x1).LT.1.0e-6)THEN
1057  intcpt = y1
1058  ELSE
1059  intcpt = y1/(slp*x1)
1060  ENDIF
1061 
1062  glb%AmpTable(1,i) = tableval(1,i)
1063  glb%AmpTable(2,i) = slp
1064  glb%AmpTable(3,i) = intcpt
1065 
1066 
1067  ENDDO
1068 
1069  DEALLOCATE(tableval)
1070 
1071  ENDIF
1072 
1073  RETURN
1074 END SUBROUTINE amplitude_sub
1075 
1076 SUBROUTINE element_sub(glb,keywd)
1077 
1078 !!****f* Rocfrac/Source/feminp/ELEMENT_SUB
1079 !!
1080 !! NAME
1081 !! ELEMENT_SUB
1082 !!
1083 !! FUNCTION
1084 !! Specifies the element type
1085 !!
1086 !! INPUTS
1087 !! glb -- global array
1088 !!
1089 !!****
1090 
1091  USE rocstar_rocfrac
1092 
1093  IMPLICIT NONE
1094 
1095  TYPE(rocfrac_global) :: glb
1096 
1097  CHARACTER(len=200) :: keywd
1098  INTEGER :: k1, k2
1099 
1100  CHARACTER(len=16) :: eltype
1101 
1102  CALL locchr(keywd,'TYPE ',4,8,k1,k2)
1103 
1104  eltype = keywd(k1:k2)
1105 
1106  SELECT CASE (trim(eltype))
1107  CASE ('V3D4')
1108  glb%iElType = 4
1109  CASE ('V3D4NCC')
1110  glb%iElType = 4
1111  glb%NdMassLump = 1
1112  glb%NdBasedEl = .true.
1113  CASE ('V3D4N')
1114  glb%iElType = 4
1115  glb%NdBasedEl = .true.
1116  CASE ('V3D10R')
1117  glb%iElType = 10
1118  glb%iElIntgratn = 1
1119  CASE ('V3D10')
1120  glb%iElType = 10
1121  CASE ('V3D10BBAR')
1122  glb%iElType = 10
1123  glb%iElIntgratn = 1
1124  CASE ('V3D8ME')
1125  glb%iElType = 8
1126  CASE default
1127  print*,' ERROR:'
1128  print*,'*ELEMENT TYPE NOT FOUND'
1129  stop
1130  END SELECT
1131 
1132  RETURN
1133 END SUBROUTINE element_sub
1134 
1135 SUBROUTINE heat_transfer_sub(glb, keywd, tmp_KappaHT, tmp_Cp, tmp_iSolnTypeHT)
1136 
1137  USE rocstar_rocfrac
1138 
1139  IMPLICIT NONE
1140 
1141  TYPE(rocfrac_global) :: glb
1142 
1143  CHARACTER :: keywd*200
1144  INTEGER :: key ! 0 = no key word found, 1 = key word found)
1145  CHARACTER*26 :: mattype
1146  INTEGER :: nummattypeht
1147  INTEGER :: i,ii,ierr, j
1148  REAL*8, DIMENSION(1:10) :: tmp_kappaht, tmp_cp
1149  INTEGER, DIMENSION(1:10) :: tmp_isolntypeht
1150 
1151  READ(glb%io_input,*) nummattypeht
1152 
1153  DO i = 1, nummattypeht
1154  glb%NumMatVolHT = glb%NumMatVolHT + 1
1155 
1156  ii = glb%NumMatVolHT
1157 
1158  IF(ii.GT.10)THEN
1159  print*,'ROCFRAC :: ERROR'
1160  print*,'Number of materials GREATER then 10'
1161  CALL mpi_finalize(glb%MPI_COMM_ROCFRAC,ierr)
1162  ENDIF
1163 
1164  READ(glb%io_input,*) j, tmp_kappaht(ii),tmp_cp(ii)
1165 
1166  tmp_isolntypeht(ii) = j ! > 0 heat transfer solution wanted, 0 = no heat transfer solution
1167 
1168 !!$! courant limit for temperature problem
1169 !!$!
1170 !!$! 2
1171 !!$! t <= Dx
1172 !!$! ------ (alpha = Kappa / (rho*cp) )
1173 !!$! 3 alpha
1174 !!$!
1175 !!$ AlphaHT = 3.* (KappaHT/RhoCp)
1176 
1177 
1178  ENDDO
1179 
1180  glb%HeatTransSoln = .true.
1181 
1182 END SUBROUTINE heat_transfer_sub
1183 
1184 
1185 
1186 SUBROUTINE probe_sub(glb)
1187 
1188 !!****f* Rocfrac/Source/feminp/PROBE_SUB
1189 !!
1190 !! NAME
1191 !! PROBE_SUB
1192 !!
1193 !! FUNCTION
1194 !! This option is used to mark a node near a
1195 !! specified coordinate
1196 !!
1197 !!
1198 !! INPUTS
1199 !! glb -- global array
1200 !!
1201 !!****
1202 
1203  USE rocstar_rocfrac
1204 
1205  IMPLICIT NONE
1206 
1207  TYPE(rocfrac_global) :: glb
1208 
1209  INTEGER :: i ! loop counter
1210 
1211  READ(glb%io_input,*) glb%NumProbesNd
1212 
1213  ALLOCATE(glb%ProbeCoorNd(1:3,1:glb%NumProbesNd))
1214  ALLOCATE(glb%ProbeNd(1:glb%NumProbesNd))
1215 
1216  DO i = 1, glb%NumProbesNd
1217  READ(glb%io_input,*) glb%ProbeCoorNd(1:3,i)
1218  ENDDO
1219 
1220  RETURN
1221 END SUBROUTINE probe_sub
1222 
1223 
1224 SUBROUTINE micromechanical_sub(glb, keywd, &
1225  tmp_e, tmp_xnu, tmp_rho, tmp_alpha, tmp_isolntype)
1226 
1227 !!****f* Rocfrac/Source/feminp/MICROMECHANICAL
1228 !!
1229 !! NAME
1230 !! MICROMECHANICAL
1231 !!
1232 !! FUNCTION
1233 !! micro mechanical material model
1234 !!
1235 !! INPUTS
1236 !! glb -- global array
1237 !! keywd -- keywd for control deck
1238 !! OUTPUT
1239 !! tmp_E -- Young's Modulus
1240 !! tmp_xnu -- Possion's ratio
1241 !! tmp_rho -- Density
1242 !! tmp_alpha -- thermal coefficient of expansion
1243 !! tmp_iSolnType -- material type model
1244 !!
1245 !!****
1246 
1247  USE rocstar_rocfrac
1248 
1249  IMPLICIT NONE
1250 
1251  TYPE(rocfrac_global) :: glb
1252 
1253  CHARACTER :: keywd*200
1254  INTEGER :: key ! 0 = no key word found, 1 = key word found)
1255  INTEGER :: k1, k2, i
1256  CHARACTER :: modeltype*16
1257  REAL*8, DIMENSION(1:10) :: tmp_e, tmp_xnu, tmp_rho, tmp_alpha
1258  INTEGER, DIMENSION(1:10) :: tmp_isolntype
1259 
1260  INTEGER :: nummatvol_loc
1261 
1262  integer :: nlgeomtype
1263 
1264  nlgeomtype = 1 ! defaults to nl
1265 
1266  CALL locchr(keywd,'MODEL ',5,8,k1,k2)
1267 
1268  IF(k1.GT.0.and.k2.GT.k1)THEN
1269 
1270  modeltype = keywd(k1:k2)
1271 
1272  IF(modeltype.EQ.'HUANG')THEN
1273  ! currently not used
1274  glb%NSTATEV = 1
1275 
1276 
1277  glb%NMATRIX = 3
1278  ALLOCATE(glb%MATRIX(1:glb%NMATRIX))
1279 
1280  READ(glb%io_input,*) glb%MATRIX(1:glb%NMATRIX), tmp_rho(1)
1281 
1282  READ(glb%io_input,*) glb%NPARTICLETYPE
1283 
1284  glb%NPARTICLE = 4
1285 
1286  allocate(glb%PARTICLE(1:glb%NPARTICLE,1:glb%NPARTICLETYPE))
1287  DO i = 1, glb%NPARTICLETYPE
1288  READ(glb%io_input,*) glb%PARTICLE(1:4,i)
1289  enddo
1290 
1291 !----- two particle sizes:
1292 ! the code in the following is limited to two particle sizes, with the
1293 ! same elastic moduli.
1294 
1295  IF (glb%particle(1,1)-glb%particle(1,2).NE.0.d0 .OR. &
1296  glb%particle(2,1)-glb%particle(2,2).NE.0.d0) THEN
1297  print*,'ROCFRAC: Error: the particles do not have the same moduli.'
1298  stop
1299  END IF
1300 
1301  IF (abs( glb%particle(3,1)+glb%particle(3,2) +glb%matrix(3)-1.d0).GT.0.001) THEN
1302  print*,'ROCFRAC: Error: the volume fractions of particles and matrix'
1303  print*,'do not add up to 1'
1304  stop
1305  END IF
1306 
1307 
1308  glb%NINTERFAC = 3
1309  allocate(glb%INTERFAC(1:glb%NINTERFAC))
1310  READ(glb%io_input,*) glb%INTERFAC(1:3)
1311 
1312  glb%NumMatVol = 1
1313  glb%DebondPart = .true.
1314 
1315 
1316  CALL homogenizedmat( tmp_rho(1), glb%Matrix(1), glb%Matrix(2), &
1317  glb%PARTICLE(1,1), glb%PARTICLE(2,1), glb%MATRIX(3), glb%cd_fastest )
1318 
1319 !!$ glb%cd_fastest = 0.d0
1320 !!$ DO i = 1, 1
1321 !!$ glb%cd_fastest = MAX( glb%cd_fastest, &
1322 !!$ SQRT(glb%MATRIX(1)*(1.d0-glb%MATRIX(2))/tmp_rho(1)/(1.d0+glb%MATRIX(2))/(1.d0-2.d0*glb%MATRIX(2))) )
1323 !!$ ENDDO
1324 !!$ DO i = 1, glb%NPARTICLETYPE
1325 !!$ glb%cd_fastest = MAX( glb%cd_fastest, &
1326 !!$ SQRT(glb%PARTICLE(1,i)*(1.d0-glb%PARTICLE(2,i))/tmp_rho(1)/(1.d0+glb%PARTICLE(2,i))/(1.d0-2.d0*glb%PARTICLE(2,i))) )
1327 !!$ ENDDO
1328 
1329 
1330  ELSE IF(modeltype.EQ.'MATOUS')THEN
1331 
1332  glb%DebondPart_Matous = .true.
1333  glb%NumMatVol = glb%NumMatVol + 1
1334 
1335  READ(glb%io_input,*) glb%NumMatVol_Part
1336 
1337  ALLOCATE(glb%E1(1:glb%NumMatVol_Part) )
1338  ALLOCATE(glb%E2(1:glb%NumMatVol_Part) )
1339  ALLOCATE(glb%E3(1:glb%NumMatVol_Part) )
1340  ALLOCATE(glb%nu12(1:glb%NumMatVol_Part) )
1341  ALLOCATE(glb%nu13(1:glb%NumMatVol_Part) )
1342  ALLOCATE(glb%nu23(1:glb%NumMatVol_Part) )
1343  ALLOCATE(glb%G12(1:glb%NumMatVol_Part) )
1344  ALLOCATE(glb%G13(1:glb%NumMatVol_Part) )
1345  ALLOCATE(glb%G23(1:glb%NumMatVol_Part) )
1346 
1347  DO i = 1, glb%NumMatVol_Part
1348  READ(glb%io_input,*) glb%E1(i), glb%E2(i), glb%E3(i), &
1349  glb%nu12(i), glb%nu13(i), glb%nu23(i), &
1350  glb%G12(i), glb%G13(i), glb%G23(i), tmp_rho(i)
1351 
1352  ENDDO
1353 
1354  READ(glb%io_input,*) glb%alpha1, glb%alpha2, glb%c2, glb%p1, glb%p2, glb%Yin, glb%a_eta, glb%a_zeta
1355 
1356  glb%cm = 1.d0 - glb%c2
1357  glb%cb = glb%c2
1358 
1359  ENDIF
1360 
1361  ELSE
1362  print*,'ERROR in *MICROMECHANICAL keywd'
1363  print*,'MODEL type not found'
1364 
1365  ENDIF
1366 
1367 
1368  tmp_isolntype(1) = nlgeomtype ! fix should not be 1
1369 
1370 
1371 END SUBROUTINE micromechanical_sub
1372 
1373 SUBROUTINE initialcondition_sub(glb,keywd)
1374 
1375 !!****f* Rocfrac/Source/feminp/INITIALCONDITION_SUB
1376 !!
1377 !! NAME
1378 !! INITIALCONDITION_SUB
1379 !!
1380 !! FUNCTION
1381 !! Specifies reference conditions
1382 !!
1383 !! INPUTS
1384 !! glb -- global array
1385 !!
1386 !!****
1387 
1388  USE rocstar_rocfrac
1389 
1390  IMPLICIT NONE
1391 
1392  TYPE(rocfrac_global) :: glb
1393 
1394  CHARACTER(len=200) :: keywd
1395  INTEGER :: k1, k2
1396 
1397  CALL locchr(keywd,'TYPE ',4,16,k1,k2)
1398 
1399 
1400  IF(keywd(k1:k2).EQ.'TEMPERATURE')THEN
1401  READ(glb%io_input,*) glb%Temperature0
1402  ENDIF
1403 
1404  RETURN
1405 END SUBROUTINE initialcondition_sub
1406 
matmodel_elastic
subroutine matmodel_elastic(glb, keywd, tmp_E, tmp_xnu, tmp_rho, tmp_alpha, tmp_iSolnType)
Definition: feminp.f90:698
element_sub
subroutine element_sub(glb, keywd)
Definition: feminp.f90:1076
dummytract_sub
subroutine dummytract_sub(glb, keywd)
Definition: feminp.f90:548
rocstar_rocfrac
Definition: ROCSTAR_RocFrac.f90:53
micromechanical_sub
subroutine micromechanical_sub(glb, keywd, tmp_E, tmp_xnu, tmp_rho, tmp_alpha, tmp_iSolnType)
Definition: feminp.f90:1224
max
Vector_n max(const Array_n_const &v1, const Array_n_const &v2)
Definition: Vector_n.h:354
conchr
subroutine conchr(text, varna, lvari, kpos0, key)
Definition: Source/CharKeyWrd.f90:149
Attribute::status
int status() const
Obtain the status of the attribute.
Definition: Attribute.h:240
matmodel_hyperelastic
subroutine matmodel_hyperelastic(glb, keywd, tmp_E, tmp_xnu, tmp_rho, tmp_alpha, tmp_iSolnType)
Definition: feminp.f90:591
sqrt
double sqrt(double d)
Definition: double.h:73
feminp
subroutine feminp(glb, myid)
Definition: feminp.f90:53
matvol_sub
subroutine matvol_sub(glb, tmp_E, tmp_xnu, tmp_rho, tmp_alpha, tmp_iSolnType)
Definition: feminp.f90:780
amplitude_sub
subroutine amplitude_sub(glb, keywd)
Definition: feminp.f90:1002
pload1_sub
subroutine pload1_sub(glb)
Definition: feminp.f90:842
Action::run
virtual void run(double t, double dt, double alpha)
Definition: Action.h:51
matcoh_sub
subroutine matcoh_sub(glb)
Definition: feminp.f90:817
i
blockLoc i
Definition: read.cpp:79
x
void int int REAL * x
Definition: read.cpp:74
initialcondition_sub
subroutine initialcondition_sub(glb, keywd)
Definition: feminp.f90:1373
boundarymm_sub
subroutine boundarymm_sub(glb)
Definition: feminp.f90:965
probe_sub
subroutine probe_sub(glb)
Definition: feminp.f90:1186
rocstar_rocfrac::rocfrac_global
Definition: ROCSTAR_RocFrac.f90:80
rchar
subroutine rchar(char, key)
Definition: Source/CharKeyWrd.f90:285
locchr
subroutine locchr(text, varna, lvari, kpos0, kpos1, kpos2)
Definition: Source/CharKeyWrd.f90:53
boundary_sub
subroutine boundary_sub(glb)
Definition: feminp.f90:888
Aff_transformation_repS2::print
virtual std::ostream & print(std::ostream &os) const
Definition: Aff_transformationS2.h:277
heat_transfer_sub
subroutine heat_transfer_sub(glb, keywd, tmp_KappaHT, tmp_Cp, tmp_iSolnTypeHT)
Definition: feminp.f90:1135
j
j indices j
Definition: Indexing.h:6
homogenizedmat
subroutine homogenizedmat(Density, Em, PRm, Ep, PRp, cm, cd_fastest)
Definition: HomogenizedMat.f90:54
prefix_sub
subroutine prefix_sub(glb)
Definition: feminp.f90:431
nodeio_sub
subroutine nodeio_sub(myid)
Definition: feminp.f90:857
key
static T_Key key
Definition: vinci_lass.c:76
dynamic_sub
subroutine dynamic_sub(glb, keywd)
Definition: feminp.f90:511
nrun_sub
subroutine nrun_sub(glb)
Definition: feminp.f90:491
boundaryht_sub
subroutine boundaryht_sub(glb)
Definition: feminp.f90:926
ale_sub
subroutine ale_sub(glb)
Definition: feminp.f90:459