Rocstar  1.0
Rocstar multiphysics simulation application
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
RocFracSubInterface.f90
Go to the documentation of this file.
1 ! * *******************************************************************
2 ! * Rocstar Simulation Suite *
3 ! * Version: Sandia Evaluation Distribution *
4 ! * Licensed To: Sandia National Laboratories *
5 ! * License Type: Evaluation *
6 ! * License Expiration: March 13, 2013 *
7 !*********************************************************************
8 ! * *******************************************************************
9 ! * Rocstar Simulation Suite *
10 ! * Copyright@2012, IllinoisRocstar LLC. All rights reserved. *
11 ! * *
12 ! * The Rocstar Simulation Suite is the property of IllinoisRocstar *
13 ! * LLC. No use or distribution of this version of the Rocstar *
14 ! * Simulation Suite beyond the license provided through separate *
15 ! * contract is permitted. *
16 ! * *
17 ! * IllinoisRocstar LLC *
18 ! * Champaign, IL *
19 ! * www.illinoisrocstar.com *
20 ! * sales@illinoisrocstar.com *
21 ! *********************************************************************
22 ! * *******************************************************************
23 ! * Initial open source Rocstar software developed by *
24 !* Center for Simulation of Advanced Rockets *
25 ! * University of Illinois at Urbana-Champaign *
26 ! * Urbana, IL *
27 !* www.csar.uiuc.edu *
28 !* *
29 ! * Copyright@2008, University of Illinois. All rights reserved. *
30 !* *
31 !* @ Redistributions of source code must retain the above copyright *
32 !* notice, this list of conditions and the following disclaimers. *
33 !* *
34 !* @ Redistributions in binary form must reproduce the above *
35 !* copyright notice, this list of conditions and the following *
36 !* disclaimers in the documentation and/or other materials *
37 !* provided with the distribution. *
38 !* *
39 !* @ Neither the names of the Center for Simulation of Advanced *
40 !* Rockets, the University of Illinois, nor the names of its *
41 !* contributors may be used to endorse or promote products derived *
42 !* from this Software without specific prior written permission. *
43 ! *********************************************************************
44 ! * *******************************************************************
45 !* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, *
46 !* EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES *
47 !* OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND *
48 !* NONINFRINGEMENT. IN NO EVENT SHALL THE CONTRIBUTORS OR *
49 !* COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER *
50 !* LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, *
51 !* ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE *
52 !* USE OR OTHER DEALINGS WITH THE SOFTWARE. *
53 !*********************************************************************
54 MODULE rocfracsubinterface
55 
56  USE rocstar_rocfracinterp
57  USE rocstar_rocfrac
58 
59  PUBLIC :: rocfracinterfaceinitial, rocfracinterfaceupdate,readsdv
60 
61 CONTAINS
62 
63  SUBROUTINE rocfracinterfaceinitial( glb,obtain_attr,surfIn)
64 
65  IMPLICIT NONE
66  include 'mpif.h'
67  include 'roccomf90.h'
68  TYPE(rocfrac_global) :: glb
69  INTEGER :: myid, numprocs, ierror
70  CHARACTER(*), INTENT(IN) :: surfin
71  INTEGER, INTENT(IN) :: obtain_attr
72  INTEGER :: npanes
73  INTEGER, POINTER, DIMENSION(:) :: paneids
74 
75  CHARACTER*9 :: timelevel
76  CHARACTER*120 :: frachdffname, meshfile
77 
78  INTEGER, POINTER :: zero, one, two
79  integer :: i,j, ierr
80  integer :: numeltypes2d
81  character(LEN=1), POINTER, DIMENSION(:) :: names
82  character(LEN=4) :: chreltype
83  integer :: endpt, startpt, chrlngth
84  integer :: pane, bcflag
85 
86  CALL mpi_comm_rank( glb%MPI_COMM_ROCFRAC, myid, ierror)
87  CALL mpi_comm_size( glb%MPI_COMM_ROCFRAC, numprocs, ierror)
88 
89 ! Subroutine to register interface data.
90  CALL mpi_barrier(glb%MPI_COMM_ROCFRAC,i)
91  IF(myid.eq.0 .AND. glb%debug_state) THEN
92  WRITE(6,'(A)') 'Rocfrac: Calling RocFrac Register...'
93  ENDIF
94 ! Get propellant density from Rocburn ! fix if more then one material
95  IF(myid.eq.0 .AND. glb%Verb.gt.1) THEN
96  WRITE(6,'(A,e11.4)') 'Rocfrac: Propellant density rho(1) is',&
97  glb%rho(1)
98  ENDIF
99 !!$ CALL COM_create_window(surWin)
100  CALL com_new_window(surwin)
101 
102 ! fix, should not this be allocated NumMaterials instead of 1 for more then one material
103  CALL com_new_attribute(surwin//'.rhos', 'w', com_double, 1, 'kg/m^3')
104  CALL com_set_array( surwin//'.rhos', 0, glb%rho,1 )
105 
106  CALL com_new_attribute(surwin//'.u', 'n', com_double, 3, 'm')
107  CALL com_new_attribute(surwin//'.vs', 'n', com_double, 3, 'm/s')
108  CALL com_new_attribute(surwin//'.uhat', 'n', com_double, 3, 'm')
109 
110 
111  CALL com_new_attribute(surwin//'.ts_alp', 'e', com_double, 1, 'Pa')
112  CALL com_new_attribute(surwin//'.bv','n',com_integer, 1, '')
113  CALL com_new_attribute( surwin//'.bf2c', 'e', com_integer, 1, '')
114  CALL com_new_attribute( surwin//'.bcflag', 'p', com_integer, 1, '')
115 
116  IF (glb%HeatTransSoln.eqv..true.) THEN
117  CALL com_new_attribute(surwin//'.qs', 'e', com_double, 1, 'W/m^2')
118  CALL com_new_attribute(surwin//'.Ts', 'n', com_double, 1, 'K')
119  ENDIF
120 
121 !!$ #OLD CALL COM_new_attribute(surWin//'.ts_alp', 'e', &
122 !!$ #OLD COM_DOUBLE_PRECISION, 3, 'Pa')
123 
124  IF ( glb%ALEenabled .eqv. .true.) CALL com_new_attribute( surwin//'.vbar_alp', 'n', com_double, 3, 'm/s')
125 
126  CALL com_get_panes( surfin, npanes, paneids)
127 
128 
129 
130  DO j = 1, npanes
131 
132  pane = paneids(j)
133 
134  CALL com_copy_array( surfin//'.bcflag', pane, bcflag)
135 
136 !!!!!!!!!!!!!!!!!!!!!! Burning Surface !!!!!!!!!!!!!!!!!!
137  IF(bcflag.EQ.1)THEN
138 
139  CALL com_get_size( surfin//".nc", pane, glb%InterfaceSFNumNodes)
140  CALL com_set_size( surwin//".nc", pane, glb%InterfaceSFNumNodes)
141  CALL com_set_size( surwin//'.bcflag', pane, 1)
142  CALL com_resize_array(surwin//'.bcflag', pane)
143 
144  CALL com_get_connectivities(surfin,pane,numeltypes2d,names)
145 
146  startpt = 1
147  DO i = 1, numeltypes2d
148  ! Search for the next attribute name
149  endpt = startpt
150  chrlngth = 0
151  DO WHILE (endpt .LE. ubound(names,1))
152  IF (names(endpt) .NE. ' ') THEN
153  chrlngth = chrlngth + 1
154  chreltype(chrlngth:chrlngth) = names(endpt)
155  endpt = endpt + 1
156  ELSE
157  EXIT
158  END IF
159  END DO
160 
161  startpt = endpt + 1
162 
163  IF(chreltype(1:chrlngth).EQ.':t6')THEN
164  CALL com_get_size( surfin//".:t6", pane, glb%InterfaceSFNumElems)
165  CALL com_set_size( surwin//".:t6", pane, glb%InterfaceSFNumElems)
166  glb%iElType2D = 6
167  ELSE IF(chreltype(1:chrlngth).EQ.':t3')THEN
168  CALL com_get_size( surfin//".:t3", pane, glb%InterfaceSFNumElems)
169  CALL com_set_size( surwin//".:t3", pane, glb%InterfaceSFNumElems)
170  glb%iElType2D = 3
171  ELSE IF(chreltype(1:chrlngth).EQ.':q4')THEN
172  CALL com_get_size( surfin//".:q4", pane, glb%InterfaceSFNumElems)
173  CALL com_set_size( surwin//".:q4", pane, glb%InterfaceSFNumElems)
174  glb%iElType2D = 4
175  ELSE
176  WRITE(0,'(A,A)') 'Rocfrac: Error: Surface mesh type',&
177  ' element not supported'
178  WRITE(0,'(A,A)') 'Read in Element Type :: ',&
179  chreltype(1:chrlngth)
180  CALL mpi_finalize(glb%MPI_COMM_ROCFRAC,ierr)
181  ENDIF
182 
183  END DO
184 
185  CALL com_free_buffer(names)
186 ! Fluids/Solids Interface Mesh
187 
188  ALLOCATE(glb%InterfaceSFNodalCoors(1:3,1:glb%InterfaceSFNumNodes))
189  ALLOCATE(glb%MapNodeSF(1:glb%InterfaceSFNumNodes))
190 
191  ALLOCATE(glb%InterfaceSFElemConn(1:glb%iElType2D,1:glb%InterfaceSFNumElems))
192  ALLOCATE(glb%MapSFElVolEl(1:glb%InterfaceSFNumElems))
193 
194 
195 ! -- The array's containing fluid-solid values
196 
197  ALLOCATE(glb%InterfaceSFNodalDisps(1:3,1:glb%InterfaceSFNumNodes))
198  ALLOCATE(glb%InterfaceSFTotalNodalDisps(1:3,1:glb%InterfaceSFNumNodes))
199  ALLOCATE(glb%InterfaceSFNodalVels(1:3,1:glb%InterfaceSFNumNodes))
200 !#OLD ALLOCATE(glb%InterfaceSFElemTract(1:3,1:glb%InterfaceSFNumElems))
201  ALLOCATE(glb%InterfaceSFElemTract(1:glb%InterfaceSFNumElems))
202 ! IF(glb%ALEenabled)THEN
203  ALLOCATE(glb%InterfaceSFVbar(1:3,1:glb%InterfaceSFNumNodes))
204 
205  IF(glb%HeatTransSoln)THEN
206  ALLOCATE(glb%InterfaceSFHeatFlux(1:glb%InterfaceSFNumElems))
207  glb%InterfaceSFHeatFlux(:) = glb%DummyFlux
208 
209  ALLOCATE(glb%InterfaceSFNodalTemp(1:glb%InterfaceSFNumNodes))
210  glb%InterfaceSFNodalTemp(:) = glb%Temperature0
211  ENDIF
212 ! ENDIF
213 
214  IF(glb%ipstatic.eqv..true.)THEN
215  ALLOCATE(glb%pstatic(1:3,1:glb%InterfaceSFNumElems))
216  glb%pstatic(1:3,1:glb%InterfaceSFNumElems) = 0.d0
217  ENDIF
218  glb%InterfaceSFNodalDisps(:,:) = 0.d0
219  glb%InterfaceSFTotalNodalDisps(:,:) = 0.d0
220  glb%InterfaceSFNodalVels(:,:) = 0.d0
221  glb%InterfaceSFElemTract(:) = 0.d0
222 !#OLD glb%InterfaceSFElemTract(:,:) = 0.d0
223  glb%InterfaceSFVbar = 0.d0
224 
225  IF ( glb%InterfaceSFNumNodes > 0) THEN ! Fluid-solid interface
226 
227 ! Register Coordinates of 2D Fluid-solid interface
228 
229 !!$ CALL COM_init_mesh( surWin//'.nc', MyId+1, glb%InterfaceSFNodalCoors, glb%InterfaceSFNumNodes)
230 
231 
232 ! print*,'glb%InterfaceSFNumNodes', glb%InterfaceSFNumNodes
233 
234  CALL com_set_array( surwin//'.nc', pane, glb%InterfaceSFNodalCoors,3 )
235 !
236 ! Registering 2D Element Connectivity of Fluid-solid interface
237 !
238  IF(glb%iElType2D.EQ.3)THEN
239 !!$ CALL COM_init_mesh( surWin//'.t3', MyId+1, glb%InterfaceSFElemConn, glb%InterfaceSFNumElems)
240 
241 
242  CALL com_set_array( surwin//'.:t3', pane, glb%InterfaceSFElemConn,3)
243 
244  ELSE IF(glb%iElType2D.EQ.6)THEN
245 !!$ CALL COM_init_mesh( surWin//'.t6', MyId+1, glb%InterfaceSFElemConn, glb%InterfaceSFNumElems)
246  CALL com_set_array( surwin//'.:t6', pane, glb%InterfaceSFElemConn,6)
247 
248  ELSE IF(glb%iElType2D.EQ.4)THEN
249 !!$ CALL COM_init_mesh( surWin//'.q4', MyId+1, glb%InterfaceSFElemConn, glb%InterfaceSFNumElems)
250 
251  CALL com_set_array( surwin//'.:q4', pane, glb%InterfaceSFElemConn,4)
252 
253  ENDIF
254 
255  CALL com_set_array( surwin//'.bv', pane, glb%MapNodeSF,1)
256  CALL com_set_array(surwin//'.bf2c', pane, glb%MapSFElVolEl, 1)
257 
258  CALL com_set_array( surwin//'.u',pane, glb%InterfaceSFNodalDisps,3)
259  CALL com_set_array( surwin//'.vs', pane, glb%InterfaceSFNodalVels,3)
260  CALL com_set_array( surwin//'.uhat', pane, glb%InterfaceSFTotalNodalDisps,3)
261 
262  CALL com_set_array( surwin//'.ts_alp', pane, glb%InterfaceSFElemTract, 1)
263 
264  IF (glb%ALEenabled) CALL com_set_array( surwin//'.vbar_alp', pane, glb%InterfaceSFVbar,3)
265 
266  IF (glb%HeatTransSoln)THEN
267  CALL com_set_array( surwin//'.qs', pane, glb%InterfaceSFHeatFlux, 1)
268  CALL com_set_array( surwin//'.Ts', pane, glb%InterfaceSFNodalTemp, 1)
269  ENDIF
270 
271  ENDIF
272 
273 !!!!!!!!!!!!!!!!!!!!!! Non-Interacting Surface !!!!!!!!!!!!!!!!!!
274  ELSE IF(bcflag.EQ.2)THEN
275 
276  CALL com_get_size( surfin//".nc", pane, glb%InterfaceSNumNodes)
277 ! WRITE(*,*) 'Number of Nodes on noninteracting structures surface: ',glb%InterfaceSNumNodes,myid
278  CALL com_set_size( surwin//".nc", pane, glb%InterfaceSNumNodes)
279 
280  CALL com_set_size( surwin//'.bcflag', pane, 1)
281  CALL com_resize_array(surwin//'.bcflag', pane)
282 
283  CALL com_get_connectivities(surfin,pane,numeltypes2d,names)
284  startpt = 1
285  DO i = 1, numeltypes2d
286  ! Search for the next attribute name
287  endpt = startpt
288  chrlngth = 0
289  DO WHILE (endpt .LE. ubound(names,1))
290  IF (names(endpt) .NE. ' ') THEN
291  chrlngth = chrlngth + 1
292  chreltype(chrlngth:chrlngth) = names(endpt)
293  endpt = endpt + 1
294  ELSE
295  EXIT
296  END IF
297  END DO
298 
299  startpt = endpt + 1
300 
301  IF(chreltype(1:chrlngth).EQ.':t6')THEN
302  CALL com_get_size( surfin//".:t6", pane, glb%InterfaceSNumElems)
303  CALL com_set_size( surwin//".:t6", pane, glb%InterfaceSNumElems)
304  glb%iElType2D = 6
305  ELSE IF(chreltype(1:chrlngth).EQ.':t3')THEN
306  CALL com_get_size( surfin//".:t3", pane, glb%InterfaceSNumElems)
307  CALL com_set_size( surwin//".:t3", pane, glb%InterfaceSNumElems)
308  glb%iElType2D = 3
309  ELSE IF(chreltype(1:chrlngth).EQ.':q4')THEN
310  CALL com_get_size( surfin//".:q4", pane, glb%InterfaceSNumElems)
311  CALL com_set_size( surwin//".:q4", pane, glb%InterfaceSNumElems)
312  glb%iElType2D = 4
313  ELSE
314  WRITE(0,'(A,A)') 'Rocfrac: Error: Surface mesh type',&
315  ' element not supported'
316  WRITE(0,'(A,A)') 'Read in Element Type :: ',&
317  chreltype(1:chrlngth)
318  CALL mpi_finalize(glb%MPI_COMM_ROCFRAC,ierr)
319  ENDIF
320 
321  END DO
322 !
323 ! Read No Solid/Fluid mesh
324 !
325  CALL com_free_buffer(names)
326 
327  ALLOCATE(glb%InterfaceSNodalCoors(1:3,1:glb%InterfaceSNumNodes))
328  ALLOCATE(glb%MapNodeS(1:glb%InterfaceSNumNodes))
329 
330  ALLOCATE(glb%InterfaceSVbar(1:3,1:glb%InterfaceSNumNodes))
331  glb%InterfaceSVbar = 0.d0
332 
333  ALLOCATE(glb%InterfaceSElemConn(1:glb%iElType2D,1:glb%InterfaceSNumElems))
334  IF(glb%EnforceTractionS.eqv..true.)THEN
335  ALLOCATE(glb%MapSElVolEl(1:glb%InterfaceSNumElems))
336  ENDIF
337 
338  IF ( glb%InterfaceSNumNodes > 0) THEN ! Non-fluid-solid interface
339 
340 ! Register Coordinates of 2D Non-Fluid-solid interface
341 
342  CALL com_set_array( surwin//'.nc', pane, glb%InterfaceSNodalCoors,3)
343  IF(glb%iElType2D.EQ.3)THEN
344 
345  CALL com_set_array( surwin//'.:t3', pane, glb%InterfaceSElemConn,3)
346 
347  ELSE IF(glb%iElType2D.EQ.6)THEN
348 
349  CALL com_set_array( surwin//'.:t6', pane, glb%InterfaceSElemConn,6)
350 
351  ELSE IF(glb%iElType2D.EQ.4)THEN
352 
353  CALL com_set_array( surwin//'.:q4', pane, glb%InterfaceSElemConn,4)
354 
355  ENDIF
356  !fix, why was this numprocs+myid+1
357 ! IF(glb%ALEenabled) CALL COM_set_array( surWin//'.vbar_alp', NumProcs+pane, glb%InterfaceSVbar,3)
358 
359 
360  IF(glb%EnforceTractionS)THEN
361  CALL com_set_array(surwin//'.bf2c', pane, glb%MapSElVolEl, 1)
362  ENDIF
363 
364 !!$ CALL COM_set_array( surWin//'.u',pane, glb%InterfaceSNodalDisps,3)
365 !!$ CALL COM_set_array( surWin//'.vs', pane, glb%InterfaceSNodalVels,3)
366 !!$ CALL COM_set_array( surWin//'.uhat', pane, glb%InterfaceSTotalNodalDisps,3)
367 
368  CALL com_set_array( surwin//'.bv', pane, glb%MapNodeS,1)
369 
370  IF( glb%ALEenabled) CALL com_set_array( surwin//'.vbar_alp', pane, glb%InterfaceSVbar,3)
371 
372  ENDIF
373 
374 
375 !!!!!!!!!!!!!!!!!!!!!! Non-Burning Interacting Surface !!!!!!!!!!!!!!!!!!
376  else if ( bcflag.eq.0 )THEN
377 
378  CALL com_get_size( surfin//".nc", pane, glb%InterfaceSFnbNumNodes)
379  CALL com_set_size( surwin//".nc", pane, glb%InterfaceSFnbNumNodes)
380 
381  CALL com_set_size( surwin//'.bcflag', pane, 1)
382  CALL com_resize_array(surwin//'.bcflag', pane)
383 
384  CALL com_get_connectivities(surfin,pane,numeltypes2d,names)
385 
386  startpt = 1
387  DO i = 1, numeltypes2d
388  ! Search for the next attribute name
389  endpt = startpt
390  chrlngth = 0
391  DO WHILE (endpt .LE. ubound(names,1))
392  IF (names(endpt) .NE. ' ') THEN
393  chrlngth = chrlngth + 1
394  chreltype(chrlngth:chrlngth) = names(endpt)
395  endpt = endpt + 1
396  ELSE
397  EXIT
398  END IF
399  END DO
400 
401  startpt = endpt + 1
402 
403  IF(chreltype(1:chrlngth).EQ.':t6')THEN
404  CALL com_get_size( surfin//".:t6", pane, glb%InterfaceSFnbNumElems)
405  CALL com_set_size( surwin//".:t6", pane, glb%InterfaceSFnbNumElems)
406  glb%iElType2D = 6
407  ELSE IF(chreltype(1:chrlngth).EQ.':t3')THEN
408  CALL com_get_size( surfin//".:t3", pane, glb%InterfaceSFnbNumElems)
409  CALL com_set_size( surwin//".:t3", pane, glb%InterfaceSFnbNumElems)
410  glb%iElType2D = 3
411  ELSE IF(chreltype(1:chrlngth).EQ.':q4')THEN
412  CALL com_get_size( surfin//".:q4", pane, glb%InterfaceSFnbNumElems)
413  CALL com_set_size( surwin//".:q4", pane, glb%InterfaceSFnbNumElems)
414  glb%iElType2D = 4
415  ELSE
416  WRITE(0,'(A,A)') 'Rocfrac: Error: Surface mesh type',&
417  ' element not supported'
418  WRITE(0,'(A,A)') 'Read in Element Type :: ',&
419  chreltype(1:chrlngth)
420  CALL mpi_finalize(glb%MPI_COMM_ROCFRAC,ierr)
421  ENDIF
422 
423  END DO
424 
425  CALL com_free_buffer(names)
426 ! Fluids/Solids Interface Mesh
427 
428  ALLOCATE(glb%InterfaceSFnbNodalCoors(1:3,1:glb%InterfaceSFnbNumNodes))
429  ALLOCATE(glb%MapNodeSFnb(1:glb%InterfaceSFnbNumNodes))
430 
431  ALLOCATE(glb%InterfaceSFnbElemConn(1:glb%iElType2D,1:glb%InterfaceSFnbNumElems))
432  ALLOCATE(glb%MapSFnbElVolEl(1:glb%InterfaceSFnbNumElems))
433 
434 
435 ! -- The array's containing fluid-solid values
436 
437  ALLOCATE(glb%InterfaceSFnbNodalDisps(1:3,1:glb%InterfaceSFnbNumNodes))
438  ALLOCATE(glb%InterfaceSFnbTotalNodalDisps(1:3,1:glb%InterfaceSFnbNumNodes))
439  ALLOCATE(glb%InterfaceSFnbNodalVels(1:3,1:glb%InterfaceSFnbNumNodes))
440 !#OLD ALLOCATE(glb%InterfaceSFElemTract(1:3,1:glb%InterfaceSFNumElems))
441  ALLOCATE(glb%InterfaceSFnbElemTract(1:glb%InterfaceSFnbNumElems))
442 ! IF(glb%ALEenabled)THEN
443  ALLOCATE(glb%InterfaceSFnbVbar(1:3,1:glb%InterfaceSFnbNumNodes))
444 ! ALLOCATE(glb%InterfaceSVbar(1:3,1:glb%InterfaceSNumNodes)) ! needed?
445 ! ENDIF
446 
447  IF(glb%ipstatic)THEN
448  ALLOCATE(glb%pstatic(1:3,1:glb%InterfaceSFNumElems))
449  glb%pstatic(1:3,1:glb%InterfaceSFNumElems) = 0.d0
450  ALLOCATE(glb%pstaticnb(1:3,1:glb%InterfaceSFnbNumElems))
451  glb%pstaticnb(1:3,1:glb%InterfaceSFnbNumElems) = 0.d0
452  ENDIF
453  glb%InterfaceSFnbNodalDisps(:,:) = 0.d0
454  glb%InterfaceSFnbTotalNodalDisps(:,:) = 0.d0
455  glb%InterfaceSFnbNodalVels(:,:) = 0.d0
456  glb%InterfaceSFnbElemTract(:) = 0.d0
457 !#OLD glb%InterfaceSFElemTract(:,:) = 0.d0
458  glb%InterfaceSFnbVbar = 0.d0
459 
460  IF ( glb%InterfaceSFnbNumNodes > 0) THEN ! Fluid-solid interface
461 
462 ! Register Coordinates of 2D Fluid-solid interface
463 
464 !!$ CALL COM_init_mesh( surWin//'.nc', MyId+1, glb%InterfaceSFNodalCoors, glb%InterfaceSFNumNodes)
465 
466 
467 ! print*,'glb%InterfaceSFNumNodes', glb%InterfaceSFNumNodes
468 
469  CALL com_set_array(surwin//'.nc', pane, glb%InterfaceSFnbNodalCoors,3 )
470 !
471 ! Registering 2D Element Connectivity of Fluid-solid interface
472 !
473  IF(glb%iElType2D.EQ.3)THEN
474 !!$ CALL COM_init_mesh( surWin//'.t3', MyId+1, glb%InterfaceSFElemConn, glb%InterfaceSFNumElems)
475 
476 
477  CALL com_set_array( surwin//'.:t3', pane, glb%InterfaceSFnbElemConn,3)
478 
479  ELSE IF(glb%iElType2D.EQ.6)THEN
480 !!$ CALL COM_init_mesh( surWin//'.t6', MyId+1, glb%InterfaceSFElemConn, glb%InterfaceSFNumElems)
481  CALL com_set_array( surwin//'.:t6', pane, glb%InterfaceSFnbElemConn,6)
482 
483  ELSE IF(glb%iElType2D.EQ.4)THEN
484 !!$ CALL COM_init_mesh( surWin//'.q4', MyId+1, glb%InterfaceSFElemConn, glb%InterfaceSFNumElems)
485 
486  CALL com_set_array( surwin//'.:q4', pane, glb%InterfaceSFnbElemConn,4)
487 
488  ENDIF
489 
490  CALL com_set_array( surwin//'.bv', pane, glb%MapNodeSFnb,1)
491  CALL com_set_array(surwin//'.bf2c', pane, glb%MapSFnbElVolEl, 1)
492 
493  CALL com_set_array( surwin//'.u',pane, glb%InterfaceSFnbNodalDisps,3)
494  CALL com_set_array( surwin//'.vs', pane, glb%InterfaceSFnbNodalVels,3)
495  CALL com_set_array( surwin//'.uhat', pane, glb%InterfaceSFnbTotalNodalDisps,3)
496 
497  CALL com_set_array( surwin//'.ts_alp', pane, glb%InterfaceSFnbElemTract, 1)
498 
499  IF ( glb%ALEenabled) CALL com_set_array( surwin//'.vbar_alp', pane, glb%InterfaceSFnbVbar,3)
500 
501  endif
502 
503  ELSE
504  WRITE(0,'(A,i4,A,i4)') 'Rocfrac: Error: Invalid bcflag', &
505  bcflag,' on surface pane',pane
506  stop
507  ENDIF
508 
509 
510  enddo
511 !!$ IF ( glb%NumNdsBC > 0) THEN ! Non-fluid-solid interface
512 !!$
513 !!$ ! CALL COM_init_mesh( surWin//'.nc', 2*NumProcs+MyId+1, glb%NumNdsBC)
514 !!$
515 !!$ CALL COM_set_size( surWin//'.nc',, MyId+1, glb%NumNdsBC)
516 !!$
517 !!$ CALL COM_set_array( surWin//'.vs', 2*NumProcs+MyId+1, glb%VeloBndry)
518 !!$ CALL COM_set_array( surWin//'.u', 2*NumProcs+MyId+1, glb%AccelBndry)
519 !!$ ENDIF
520 
521 ! No longer needed
522 !!$ CALL COM_new_attribute( surWin//'.bcflag', 'p', COM_INTEGER, 1, '')
523 !!$ CALL COM_allocate_array( surWin//'.bcflag', MyId+1)
524 !!$
525 !!$ stop
526 !!$
527 !!$ IF ( glb%InterfaceSFNumNodes > 0) THEN ! Fluid-solid interface
528 !!$ CALL COM_get_array( surWin//'.bcflag', MyId+1, one)
529 !!$ one = 1
530 !!$ END IF
531 !!$
532 !!$ IF ( glb%InterfaceSNumNodes > 0) THEN ! Non-fluid-solid interface
533 !!$ CALL COM_get_array( surWin//'.bcflag', NumProcs+MyId+1, two)
534 !!$ two = 2
535 !!$ ENDIF
536 ! put into volume mesh
537 
538 !!$ IF ( glb%NumNdsBC > 0) THEN ! Surface Boundary Condition flags
539 !!$ CALL COM_get_array( surWin//'.bcflag', 2*NumProcs+MyId+1, two)
540 !!$ two = 2
541 !!$ ENDIF
542 
543 
544  CALL com_window_init_done( surwin)
545 
546 
547 
548  CALL com_free_buffer(paneids)
549 
550  CALL com_call_function( obtain_attr, 2, &
551  com_get_attribute_handle_const( surfin//".all"), &
552  com_get_attribute_handle( surwin//".all"))
553 
554  CALL mpi_barrier(glb%MPI_COMM_ROCFRAC,i)
555  IF(myid.eq.0 .AND. glb%debug_state) THEN
556  WRITE(6,'(A)') 'Rocfrac: Calling RocFrac Register... Done'
557  ENDIF
558 11 FORMAT(a,'_',a,a1)
559 !---------------------------------------------------------------------------------------------------
560 
561  END SUBROUTINE rocfracinterfaceinitial
562 
563 
564  SUBROUTINE rocfracinterfaceupdate( glb, CurrentTime, &
565  currenttimestep, man_update_inbuff)
566  IMPLICIT NONE
567 
568  TYPE(rocfrac_global), POINTER :: glb
569  REAL*8, INTENT(IN) :: currenttime, currenttimestep
570  INTEGER, INTENT(IN) :: man_update_inbuff
571 
572 ! Subroutine to gather current interface position computed by the solids code into
573 ! the database used by the interpolation procedures.
574 
575  ! Local variables:
576 
577  REAL*8 :: solidsdisp, solidscoor
578 
579  INTEGER :: iinterfacenode, solidsnodenum
580 
581 ! Update ALL the surface mesh coordinates
582 
583 ! Burning/Interacting surfaces
584 
585  DO iinterfacenode = 1, glb%InterfaceSFNumNodes
586 
587  solidsnodenum = abs(glb%MapNodeSF(iinterfacenode))
588 
589  glb%InterfaceSFNodalCoors(1,iinterfacenode) = glb%MeshCoor(1,solidsnodenum )
590  glb%InterfaceSFNodalCoors(2,iinterfacenode) = glb%MeshCoor(2,solidsnodenum )
591  glb%InterfaceSFNodalCoors(3,iinterfacenode) = glb%MeshCoor(3,solidsnodenum )
592 
593  END DO
594 
595 ! Non-Burning/Interacting surfaces
596 
597  DO iinterfacenode = 1, glb%InterfaceSFnbNumNodes
598 
599  solidsnodenum = abs(glb%MapNodeSFnb(iinterfacenode))
600 
601  glb%InterfaceSFnbNodalCoors(1,iinterfacenode) = glb%MeshCoor(1,solidsnodenum )
602  glb%InterfaceSFnbNodalCoors(2,iinterfacenode) = glb%MeshCoor(2,solidsnodenum )
603  glb%InterfaceSFnbNodalCoors(3,iinterfacenode) = glb%MeshCoor(3,solidsnodenum )
604 
605  END DO
606 
607 ! Non-Interacting surfaces
608 
609  DO iinterfacenode = 1, glb%InterfaceSNumNodes
610 
611  solidsnodenum = abs(glb%MapNodeS(iinterfacenode))
612  IF(solidsnodenum .ne. 0) THEN
613  glb%InterfaceSNodalCoors(1,iinterfacenode) = glb%MeshCoor(1,solidsnodenum )
614  glb%InterfaceSNodalCoors(2,iinterfacenode) = glb%MeshCoor(2,solidsnodenum )
615  glb%InterfaceSNodalCoors(3,iinterfacenode) = glb%MeshCoor(3,solidsnodenum )
616  ENDIF
617 
618  END DO
619 !---------------------------------------------------------------------------------------------------
620 
621  END SUBROUTINE rocfracinterfaceupdate
622 
623 
624  SUBROUTINE readsdv(glb,myid)
625 
626  IMPLICIT NONE
627 
628  include 'roccomf90.h'
629  include 'mpif.h'
630 
631  INTEGER :: myid
632 
633  TYPE(rocfrac_global),POINTER :: glb
634 
635  INTEGER :: pid
636 
637 
638  INTEGER :: hdl_read, hdl_obtain, hdl_all
639 
640  CHARACTER(*), PARAMETER :: overlaywin = "Overlay"
641  CHARACTER(*), PARAMETER :: overlaywin2 = "Overlay2"
642 
643 
644  INTEGER :: comm_self
645  CHARACTER(*), PARAMETER :: prefix1 = "A"
646  CHARACTER(*), PARAMETER :: prefix2 = "B"
647 
648  CHARACTER(LEN=5) :: sdv_material
649  CHARACTER(LEN=12) :: sdv_wname
650  CHARACTER(LEN=27) :: fname1, fname2
651  CHARACTER(LEN=3) :: ichr3
652  INTEGER :: i,j
653 
654 
655  INTEGER, POINTER, DIMENSION(:) :: mapfaceel2vol1a, faceofvolel1a
656 
657  INTEGER :: iprocs,ios, iaux
658 ! obtain function handle ------------------------------------------------------
659 
660 
661 ! IF(myid.eq.1) THEN
662 
663  pid = (myid+1)*100 + 3
664 
665 
666  WRITE(ichr3,'(I3.3)') pid ! problem if over 999 processors, fix
667 
668  fname1 = 'Rocfrac/Rocin/A_'//ichr3//'_sdv.hdf'
669  fname2 = 'Rocfrac/Rocin/B_'//ichr3//'_sdv.hdf'
670 
671 
672  CALL rocin_load_module( "SDV_IN")
673 
674  hdl_read = com_get_function_handle( 'SDV_IN.read_windows')
675  hdl_obtain = com_get_function_handle( 'SDV_IN.obtain_attribute')
676 
677 
678 
679 
680 ! Define the base-window and sdv-window names
681 
682  sdv_material = prefix1//'_sdv'
683  sdv_wname = overlaywin//sdv_material
684  CALL com_new_window( overlaywin )
685 
686 ! // Read the pane from the given file. Note that the file contains both
687 ! // the parent and the subdivided windows. Read only the subdivided one.
688  comm_self = mpi_comm_self
689  CALL com_call_function( hdl_read, 4, fname1, overlaywin, &
690  sdv_material, comm_self)
691  hdl_all = com_get_attribute_handle( sdv_wname//".all")
692  CALL com_call_function( hdl_obtain, 3, hdl_all, hdl_all, pid)
693 ! // Obtain number of sub-nodes, sub-faces, nodes, and faces
694 
695  CALL com_get_size( sdv_wname//'.sn_parent_fcID', pid, glb%nsubn1)
696 
697 ! PRINT*,'Number of sub-nodes =',glb%nsubn1
698 
699  CALL com_get_size( sdv_wname//'.:t3', pid, glb%nsubf1)
700 
701 ! PRINT*,'Number of sub-faces =', glb%nsubf1
702 
703  ALLOCATE(glb%Sthresh1(1:3, 1:glb%nsubf1))
704  glb%Sthresh1(1:3, 1:glb%nsubf1) = glb%Sinit(1) ! fix, this should be read-in
705 
706 
707 ! // Obtain the connectivity
708 
709  ALLOCATE( glb%sd_subfaces1(1:3,1:glb%nsubf1) )
710 
711  CALL com_copy_array( sdv_wname//".:t3", pid, glb%sd_subfaces1, 3)
712 
713 !!$ DO i = 1, glb%nsubf1
714 !!$ PRINT*,glb%sd_subfaces1(1:3,i)
715 !!$ ENDDO
716 
717  ALLOCATE(glb%sd_coor1(1:3,1:glb%nsubn1))
718 
719  CALL com_copy_array( sdv_wname//".nc", pid, glb%sd_coor1, 3)
720 
721 !!$ DO i = 1, glb%nsubn1
722 !!$ PRINT*,glb%sd_coor1(1:3,i)
723 !!$ ENDDO
724 
725  ALLOCATE(glb%sd_subface_parents1(1:glb%nsubf1))
726 
727  CALL com_copy_array( sdv_wname//".sf_parent", pid, glb%sd_subface_parents1)
728 
729 !!$ DO i = 1, glb%nsubf1
730 !!$ PRINT*,'sd_subface_parents',glb%sd_subface_parents1(i)
731 !!$ ENDDO
732 
733 
734  ALLOCATE(glb%sd_subface_nat_coors1(1:6,1:glb%nsubf1))
735 
736 
737 
738  CALL com_copy_array( sdv_wname//".sf_ncs", pid, &
739  glb%sd_subface_nat_coors1, 6)
740 
741 !
742 ! // NOTE: The last argument (6) indicates that the local coordinates are
743 ! // stored consecutively (xi1, eta1, xi2, eta2, xi3, eta3). Use the number
744 ! // one (1) if the xi1 for all nodes are stored together and then xi2, etc..
745 
746 !!$ DO i = 1, glb%nsubf1
747 !!$ PRINT*, i,glb%nsubf1,glb%sd_subface_nat_coors1(1:6,i)
748 !!$ END DO
749  ALLOCATE(glb%sd_subface_counterparts1(1:glb%nsubf1) )
750 
751 
752  CALL com_copy_array( sdv_wname//".sf_cntpt_fcID", pid, glb%sd_subface_counterparts1)
753 
754 
755 
756 !!$ DO i = 1, glb%nsubf1
757 !!$ PRINT*,glb%sd_subface_counterparts1(i)
758 !!$ ENDDO
759 
760 
761 
762 !!$ PRINT*,'finished readsdv',myid
763 !!$ CALL MPI_BARRIER(glb%MPI_COMM_ROCFRAC,i)
764 !!$ stop
765 
766 ! // Delete the window created by Rocin.
767 
768 !!$ CALL COM_delete_window( sdv_wname)
769 !!$
770 !!$
771 !!$ CALL Rocin_unload_module( "SDV_IN")
772 !!$
773 !!$
774 !!$ CALL Rocin_load_module( "SDV_IN")
775 
776 ! Define the base-window and sdv-window names
777 
778  sdv_material = prefix2//'_sdv'
779  sdv_wname = overlaywin//sdv_material
780 
781 ! CALL COM_new_window( OverlayWin )
782 
783 
784 
785 
786 
787 ! // Read the pane from the given file. Note that the file contains both
788 ! // the parent and the subdivided windows. Read only the subdivided one.
789  comm_self = mpi_comm_self
790  CALL com_call_function( hdl_read, 4, fname2, overlaywin, &
791  sdv_material, comm_self)
792  hdl_all = com_get_attribute_handle( sdv_wname//".all")
793  CALL com_call_function( hdl_obtain, 3, hdl_all, hdl_all, pid)
794  CALL com_get_size( sdv_wname//'.sn_parent_fcID', pid, glb%nsubn2)
795 
796 
797 
798 
799 ! PRINT*,'Number of sub-nodes =',glb%nsubn2
800 
801  CALL com_get_size( sdv_wname//'.:t3', pid, glb%nsubf2)
802 
803 
804 
805 ! PRINT*,'Number of sub-faces =', glb%nsubf2
806 
807  ALLOCATE(glb%sd_subface_counterparts2(1:glb%nsubf2) )
808 
809  CALL com_copy_array( sdv_wname//".sf_cntpt_fcID", pid, glb%sd_subface_counterparts2)
810 
811 
812 ! // Obtain the connectivity
813 
814  ALLOCATE( glb%sd_subfaces2(1:3,1:glb%nsubf1) ) ! GETS STUCK HERE FOR more then ONE processor
815 
816 !!$ PRINT*,'finished readsdv!!!!!!!!!',myid,glb%nsubn2,glb%nsubf2
817 !!$ CALL MPI_BARRIER(glb%MPI_COMM_ROCFRAC,i)
818 
819  CALL com_copy_array( sdv_wname//".:t3", pid, glb%sd_subfaces2, 3)
820 
821 !!$ DO i = 1, glb%nsubf2
822 !!$ PRINT*,glb%sd_subfaces2(1:3,i)
823 !!$ ENDDO
824 
825  ALLOCATE(glb%sd_coor2(1:3,1:glb%nsubn2))
826  CALL com_copy_array( sdv_wname//".nc", pid, glb%sd_coor2, 3)
827 
828 
829 !!$ DO i = 1, glb%nsubn2
830 !!$ PRINT*,glb%sd_coor2(1:3,i)
831 !!$ ENDDO
832 
833 
834 
835 ! PRINT*,'ldfdfkjsdflkj'
836  ALLOCATE(glb%sd_subface_nat_coors2(1:6,1:glb%nsubf2))
837 ! PRINT*,'dsflkjsdflkj'
838  CALL com_copy_array( sdv_wname//".sf_ncs", pid, &
839  glb%sd_subface_nat_coors2, 6)
840 !
841 ! PRINT*,'lkjsdflkj'
842 
843 !
844 ! // NOTE: The last argument (6) indicates that the local coordinates are
845 ! // stored consecutively (xi1, eta1, xi2, eta2, xi3, eta3). Use the number
846 ! // one (1) if the xi1 for all nodes are stored together and then xi2, etc..
847 
848 !!$ DO i = 1, glb%nsubf2
849 !!$ PRINT*, i,glb%nsubf2,glb%sd_subface_nat_coors2(1:6,i)
850 !!$ END DO
851 
852 
853 
854 !!$ DO i = 1, glb%nsubf2
855 !!$ PRINT*,i,glb%nsubf2,glb%sd_subface_counterparts2(i)
856 !!$ ENDDO
857 
858 ! PRINT*,'lkjsdflkj',glb%nsubf2,myid
859  ALLOCATE(glb%sd_subface_parents2(1:glb%nsubf2))
860 ! PRINT*,'copy',myid
861  CALL com_copy_array( sdv_wname//".sf_parent", pid, glb%sd_subface_parents2)
862 ! PRINT*,'sdfsdflkjsdflkj',myid
863 !!$ DO i = 1, glb%nsubf2
864 !!$ PRINT*,'sd_subface_parents',glb%sd_subface_parents2(i)
865 !!$ ENDDO
866 
867 
868 ! // Delete the window created by Rocin.
869  CALL com_delete_window( sdv_wname)
870 ! PRINT*,'sdfsdflkjsdflkj',myid
871 
872 ! // Unload Rocin from Roccom.
873 
874  CALL rocin_unload_module( "SDV_IN")
875 
876 
877 
878 ! CALL COM_new_attribute( surWin//'.bf2c', 'e', COM_INTEGER, 1, '')
879 ! CALL COM_allocate_array(surWin//'.bf2c', iNI, ElFlag_List, 1)
880 ! CALL COM_new_attribute( surWin//'.faceOnCell', 'e', COM_INTEGER, 1, '')
881 ! CALL COM_allocate_array(surWin//'.faceOnCell', iNI, FaceOnCell, 1)
882 
883 
884 ! ENDIF
885 
886 ! stop
887 
888 
889 
890 END SUBROUTINE readsdv
891 
892 END MODULE rocfracsubinterface
893 
Mesquite::Matrix3D::zero
void zero()
Sets all entries to zero (more efficient than assignement).
Definition: includeLinks/Matrix3D.hpp:167
rocstar_rocfrac
Definition: ROCSTAR_RocFrac.f90:53
rocfracsubinterface::rocfracinterfaceupdate
subroutine, public rocfracinterfaceupdate(glb, CurrentTime, CurrentTimeStep, MAN_update_inbuff)
Definition: RocFracSubInterface.f90:564
rocstar_rocfracinterp
Definition: ROCSTAR_RocFracInterp.f90:53
i
blockLoc i
Definition: read.cpp:79
rocstar_rocfrac::rocfrac_global
Definition: ROCSTAR_RocFrac.f90:80
Attribute::pane
const Pane * pane() const
Obtain a constant pointer to the owner pane of the attribute.
Definition: Attribute.h:172
j
j indices j
Definition: Indexing.h:6
rocfracsubinterface::rocfracinterfaceinitial
subroutine, public rocfracinterfaceinitial(glb, obtain_attr, surfIn)
Definition: RocFracSubInterface.f90:63
readsdv
program readsdv
Definition: readsdv.f90:23
Line_2::a
RT a() const
Definition: Line_2.h:140