Rocstar  1.0
Rocstar multiphysics simulation application
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
utilities/RocfracPrep/2Dmesh.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 SUBROUTINE mesh2d(nprocs,iProcs,ichr4)
55 
56  use meshdata
57  USE linked_list2
58 
59  IMPLICIT NONE
60 
61  include 'roccomf90.h'
62 
63  CHARACTER*4 :: ichr4
64  INTEGER :: new
65  INTEGER :: i,k
66  INTEGER :: icounter
67  INTEGER :: iprocs
68  INTEGER :: nprocs
69  INTEGER :: numnpnew
70  INTEGER :: glbndnum
71  LOGICAL :: surfaceelonproc
72  INTEGER :: testnumsurfel
73  INTEGER, POINTER, DIMENSION(:,:) :: nodeflag2d
74  INTEGER :: itestcnt
75 
76 ! User-defined list element
77 ! The Link_Type field MUST be the FIRST in the user-defined list element
78 ! Note pointer to data so as to easily create sublists
79 
80  TYPE user_type_surfndlist
81  TYPE(link_type) :: link
82  TYPE(user_data_type_procnodelist), POINTER :: data
83  END TYPE user_type_surfndlist
84 
85  TYPE user_data_type_surfndlist
86  INTEGER :: glbndnum
87  END TYPE user_data_type_surfndlist
88 
89 ! Auxilliary data type required for the transfer function
90  TYPE user_ptr_type_surfndlist
91  TYPE(user_type_surfndlist), POINTER :: p
92  END TYPE user_ptr_type_surfndlist
93 
94  TYPE(user_ptr_type_surfndlist) :: user_surfndlist
95  TYPE(list_type) :: surfndlist
96 
97  TYPE(link_ptr_type) :: link
98 
99  TYPE(surfmesh_tri3_ptr),POINTER :: ptrtri3
100  TYPE(surfmesh_tri6_ptr), POINTER :: ptrtri6
101  TYPE(surfmesh_hex8_ptr), POINTER :: ptrhex8
102 
103  INTEGER, POINTER, DIMENSION(:) :: bcsurfflagzero,bcsurfflagone,bcsurfflagtwo
104 
105 
106  CHARACTER(*), PARAMETER :: surwin = "sfrac"
107 
108  REAL*8, POINTER, DIMENSION(:,:) :: meshcoor
109 
110  INTEGER, POINTER, DIMENSION(:) :: elflag_list
111 
112  INTEGER :: write_attr, set_option, ierrflg, sur_all
113  INTEGER :: inb, ibu, ini
114  ! MS : A toggle to check whether append or write
115  LOGICAL :: appendtoggle = .false.
116  appendtoggle = .false.
117  WRITE (*,*) '**** iProc = ', iprocs, ' appendToggle = ', appendtoggle
118  ! MS
119 
120 ! obtain function handle ------------------------------------------------------
121 
122  write_attr = com_get_function_handle( 'OUT.write_attribute')
123  set_option = com_get_function_handle( 'OUT.set_option')
124 
125 ! do not append process rank -----------------
126 
127  CALL com_call_function( set_option, 2, 'rankwidth', '0')
128 ! CALL COM_call_function( set_option, 2, 'mode', 'w')
129 
130  CALL com_new_window( surwin )
131 
132 ! --------------------------------------------------
133 ! ------------------------------------------------
134 ! (1) Ignitable interface (BURNING)
135 ! ---------------------------------------
136 ! -------------------------------------
137 
138  inb = iprocs*100 + 1
139  ibu = iprocs*100 + 2
140  ini = iprocs*100 + 3
141 
142  numnpnew = 0
143 
144  IF(meshtype2d.EQ.6)THEN
145 
146 ! renumber elements nodes
147 
148  ptrtri6 => surfmesh_tri6_sf_head
149  DO WHILE(ASSOCIATED(ptrtri6).eqv..true.)
150 ! mark that changed surface mesh by changing it to a negative
151  DO i = 1, 6
152  glbndnum = ptrtri6%ElemData(i)
153  IF( nodeflag(glbndnum).GT.0) THEN ! if zero means not on this processor, neg means already renumber node
154  numnpnew = numnpnew + 1
155  nodeflag(glbndnum) = -nodeflag(glbndnum) ! overwrite nodeflag here but i need it for non S/F interface
156  ENDIF
157  ENDDO
158  ptrtri6 => ptrtri6%next
159  ENDDO
160 
161 ! check for special situation where there is just an edge on the surface
162 
163  ptrtri6 => surfmesh_tri6_sf_head
164  DO WHILE(ASSOCIATED(ptrtri6).eqv..true.)
165  DO i = 1, 6
166  glbndnum = ptrtri6%ElemData(i)
167  DO k = 1,maxnumberofprocstosharenode
168  IF(procndlist(glbndnum,k).EQ.iprocs.AND.nodeflag(glbndnum).GT.0)THEN
169  print*,'Detected surface sliver for pressure'
170  numnpnew = numnpnew + 1
171  nodeflag(glbndnum) = -nodeflag(glbndnum) ! overwrite nodeflag here but i need it for non S/F interface
172  ENDIF
173  ENDDO
174  ENDDO
175  ptrtri6 => ptrtri6%next
176  ENDDO
177 
178 ! write the surface mesh files
179 
180  IF(numeltet2d(1,iprocs).NE.0.AND.numnpnew.EQ.0)THEN
181  print*,'ERROR, found in surface mesh for Solid/Fluid Ignitable interface'
182  print*,'Found ',numeltet2d(1,iprocs),'triangles but no nodes'
183  print*,'Stopping'
184  stop
185  ENDIF
186 
187  IF(numeltet2d(1,iprocs).EQ.0.AND.numnpnew.NE.0)THEN
188  print*,'Warning: found in sufrace mesh for Solid/Fluid Ignitable interface'
189  print*,'Found ',numnpnew,' nodes but no triangles'
190  ENDIF
191 
192 
193 ! IF(NumNPNew.NE.0.AND.NumEltet2D(1,iProcs).NE.0) &
194 ! write(4003,*) 'ZONE N=', NumNPNew, 'E=',NumEltet2D(1,iProcs)*4,'F=FEPOINT, ET=TRIANGLE'
195 
196  nodeflag(:)= abs( nodeflag(:)) ! reset
197 
198  ALLOCATE(nodeflag2d(1:2,1:numnpnew)) ! for storing what I'm overwriting
199 
200  numnpnew = 0
201 ! renumber elements nodes
202 
203  ptrtri6 => surfmesh_tri6_sf_head
204  DO WHILE(ASSOCIATED(ptrtri6).eqv..true.)
205 ! mark that changed surface mesh by changing it to a negative
206  DO i = 1, 6
207  glbndnum = ptrtri6%ElemData(i)
208  IF( nodeflag(glbndnum).GT.0) THEN ! if zero means not on this processor
209  numnpnew = numnpnew + 1
210 !!$ WRITE(4001,'(3(x,e16.9),1i10)') coor(1:3,glbNdNum),NodeFlag(glbNdNum)
211 !!$ WRITE(4003,*) coor(1:3,glbNdNum)
212  nodeflag2d(1,numnpnew) = glbndnum
213  nodeflag2d(2,numnpnew) = nodeflag(glbndnum)
214  nodeflag(glbndnum) = - numnpnew
215  ENDIF
216  ENDDO
217  ptrtri6 => ptrtri6%next
218  ENDDO
219  ptrtri6 => surfmesh_tri6_sf_head
220  DO WHILE(ASSOCIATED(ptrtri6).eqv..true.)
221  DO i = 1, 6
222  glbndnum = ptrtri6%ElemData(i)
223  DO k = 1,maxnumberofprocstosharenode
224  IF(procndlist(glbndnum,k).EQ.iprocs.AND.nodeflag(glbndnum).GT.0)THEN
225  numnpnew = numnpnew + 1
226  print*,'ROCFRAC: Dectected surface sliver for pressure, writing'
227 ! WRITE(4001,'(3(x,e16.9),1i10)') coor(1:3,glbNdNum),NodeFlag(glbNdNum)
228  nodeflag2d(1,numnpnew) = glbndnum
229  nodeflag2d(2,numnpnew) = nodeflag(glbndnum)
230  nodeflag(glbndnum) = - numnpnew
231  ENDIF
232  ENDDO
233  ENDDO
234  ptrtri6 => ptrtri6%next
235  ENDDO
236 
237  CALL com_new_attribute( surwin//'.nc', 'n', com_double, 3, 'm')
238  CALL com_set_size( surwin//'.nc', ibu, numnpnew )
239  CALL com_resize_array(surwin//'.nc', ibu, meshcoor, 3)
240 
241  DO i = 1, numnpnew
242  meshcoor(1:3,i) = coor(1:3,nodeflag2d(1,i))
243  END DO
244 
245 ! nodeflag to VolumeNode
246 
247  CALL com_new_attribute( surwin//'.bv', 'n', com_integer, 1, 'm')
248  CALL com_set_array( surwin//'.bv', ibu, nodeflag2d(2,1), 2)
249 
250 
251 ! connectiveity
252  testnumsurfel = numeltet2d(1,iprocs)
253 
254  CALL com_set_size( surwin//'.:t6', ibu, testnumsurfel)
255  CALL com_resize_array( surwin//'.:t6', ibu, elconntable, 6)
256 
257 
258 ! elemflag to VolumeElem
259 
260  CALL com_new_attribute( surwin//'.bf2c', 'e', com_integer, 1, '')
261  CALL com_resize_array(surwin//'.bf2c', ibu, elflag_list, 1)
262 
263  CALL com_new_attribute( surwin//'.bcflag', 'p', com_integer, 1, '')
264  CALL com_set_size( surwin//'.bcflag', ibu, 1)
265  CALL com_resize_array(surwin//'.bcflag', ibu, bcsurfflagone, 1)
266 
267  bcsurfflagone(1) = 1
268 
269  testnumsurfel = 0
270  ptrtri6 => surfmesh_tri6_sf_head
271  DO WHILE(ASSOCIATED(ptrtri6).eqv..true.)
272  surfaceelonproc = .false.
273 !!$ do i = 1, 6
274 !!$ IF(NodeFlag(ptrtri6%ElemData(i)).eq.0) THEN
275 !!$ SurfaceElOnProc = .false.
276 !!$ exit
277 !!$ endif
278 !!$ enddo
279 
280  IF(epart(ptrtri6%ElemData(7)).EQ.iprocs) surfaceelonproc = .true.
281  IF(surfaceelonproc) THEN
282 !!$ WRITE(4001,'(10i10)') ABS(NodeFlag(ptrtri6%ElemData(1:6))),&
283 !!$ ElFlag(ptrtri6%ElemData(7)),1,1 ! last two are not really used
284 
285  testnumsurfel = testnumsurfel + 1
286  elconntable(1:6,testnumsurfel) = abs(nodeflag(ptrtri6%ElemData(1:6)))
287  elflag_list(testnumsurfel) = elflag(ptrtri6%ElemData(7))
288 
289 !!$ write(4003,*) ABS(NodeFlag(ptrtri6%ElemData(1))),ABS(NodeFlag(ptrtri6%ElemData(4))),ABS(NodeFlag(ptrtri6%ElemData(6)))
290 !!$ write(4003,*) ABS(NodeFlag(ptrtri6%ElemData(4))),ABS(NodeFlag(ptrtri6%ElemData(2))),ABS(NodeFlag(ptrtri6%ElemData(5)))
291 !!$ write(4003,*) ABS(NodeFlag(ptrtri6%ElemData(6))),ABS(NodeFlag(ptrtri6%ElemData(4))),ABS(NodeFlag(ptrtri6%ElemData(5)))
292 !!$ write(4003,*) ABS(NodeFlag(ptrtri6%ElemData(3))),ABS(NodeFlag(ptrtri6%ElemData(6))),ABS(NodeFlag(ptrtri6%ElemData(5)))
293  ENDIF
294  ptrtri6 => ptrtri6%next
295  ENDDO
296 
297  IF(testnumsurfel.NE.numeltet2d(1,iprocs))THEN
298  print*,'ERROR, number of triangles from link list in for ignitable surfaces'
299  print*,' different then that in read_patran'
300  print*,testnumsurfel,numeltet2d(1,iprocs)
301  stop
302  ENDIF
303 
304  nullify(ptrtri6)
305 
306 
307  ! obtain function handle ------------------------------------------------------
308 
309 ! write_attr = COM_get_function_handle( 'OUT.write_attribute')
310 ! set_option = COM_get_function_handle( 'OUT.set_option')
311 
312 ! IF(iProcs.EQ.1)THEN
313 
314  CALL com_call_function( set_option, 2, 'mode', 'w')
315 
316 ! ELSE
317 !
318 ! CALL COM_call_function( set_option, 2, 'mode', 'a')
319 !
320 ! ENDIF
321 
322 ! do not append process rank -----------------
323 
324 ! CALL COM_call_function( set_option, 2, 'rankwidth', '0')
325 
326 ! write surface window ------------------------
327 
328 ! sur_all = Com_get_attribute_handle( surWin//'.all')
329 ! CALL COM_call_function( write_attr, 4, 'Rocin/SurfMesh_FS_Ignt.hdf', sur_all,&
330 ! "solid_surf","00.000000")
331 
332 ! CALL COM_delete_window( surWin)
333 
334 
335  ELSE IF(meshtype2d.EQ.3)THEN
336 
337 ! renumber elements nodes
338 
339  ptrtri3 => surfmesh_tri3_sf_head
340  do while(associated(ptrtri3))
341 ! mark that changed surface mesh by changing it to a negative
342  DO i = 1, 3
343  glbndnum = ptrtri3%ElemData(i)
344  IF( nodeflag(glbndnum).GT.0) THEN ! if zero means not on this processor, neg means already renumber node
345  numnpnew = numnpnew + 1
346  nodeflag(glbndnum) = -nodeflag(glbndnum) ! overwrite nodeflag here but i need it for non S/F interface
347  endif
348  enddo
349  ptrtri3 => ptrtri3%next
350  ENDDO
351 
352 ! check for special situation where there is just an edge on the surface
353 
354  ptrtri3 => surfmesh_tri3_sf_head
355  do while(associated(ptrtri3))
356  DO i = 1, 3
357  glbndnum = ptrtri3%ElemData(i)
358  DO k = 1,maxnumberofprocstosharenode
359  IF(procndlist(glbndnum,k).EQ.iprocs.AND.nodeflag(glbndnum).GT.0)THEN
360  print*,'Dectected surface sliver for pressure'
361  numnpnew = numnpnew + 1
362  nodeflag(glbndnum) = -nodeflag(glbndnum) ! overwrite nodeflag here but i need it for non S/F interface
363  ENDIF
364  ENDDO
365  ENDDO
366  ptrtri3 => ptrtri3%next
367  ENDDO
368 
369 ! write the surface mesh files
370 
371  IF(numeltet2d(1,iprocs).NE.0.AND.numnpnew.EQ.0)THEN
372  print*,'ERROR, found in surface mesh for Solid/Fluid interface'
373  print*,'Found ',numeltet2d(1,iprocs),'triangles but no nodes'
374  print*,'Stopping'
375  stop
376  ENDIF
377 
378  IF(numeltet2d(1,iprocs).EQ.0.AND.numnpnew.NE.0)THEN
379  print*,'Warning: found in sufrace mesh for Solid/Fluid interface'
380  print*,'Found ',numnpnew,' nodes but no triangles'
381  ENDIF
382 
383  nodeflag(:)= abs( nodeflag(:)) ! reset
384 
385  ALLOCATE(nodeflag2d(1:2,1:numnpnew)) ! for storing what I'm overwriting
386 
387  numnpnew = 0
388 ! renumber elements nodes
389 
390  ptrtri3 => surfmesh_tri3_sf_head
391  do while(associated(ptrtri3).eqv..true.)
392 ! mark that changed surface mesh by changing it to a negative
393  DO i = 1, 3
394  glbndnum = ptrtri3%ElemData(i)
395  IF( nodeflag(glbndnum).GT.0) THEN ! if zero means not on this processor
396  numnpnew = numnpnew + 1
397 !!$ WRITE(4001,'(3(x,e16.9),1i10)') coor(1:3,glbNdNum),NodeFlag(glbNdNum)
398 !!$ write(4003,*) coor(1:3,glbNdNum)
399  nodeflag2d(1,numnpnew) = glbndnum
400  nodeflag2d(2,numnpnew) = nodeflag(glbndnum)
401  nodeflag(glbndnum) = - numnpnew
402  endif
403  enddo
404  ptrtri3 => ptrtri3%next
405  ENDDO
406  ptrtri3 => surfmesh_tri3_sf_head
407  do while(associated(ptrtri3).eqv..true.)
408  DO i = 1, 3
409  glbndnum = ptrtri3%ElemData(i)
410  DO k = 1,maxnumberofprocstosharenode
411  IF(procndlist(glbndnum,k).EQ.iprocs.AND.nodeflag(glbndnum).GT.0)THEN
412  numnpnew = numnpnew + 1
413  print*,'Dectected surface sliver for pressure, writing'
414 ! WRITE(4001,'(3(x,e16.9),1i10)') coor(1:3,glbNdNum),NodeFlag(glbNdNum)
415  nodeflag2d(1,numnpnew) = glbndnum
416  nodeflag2d(2,numnpnew) = nodeflag(glbndnum)
417  nodeflag(glbndnum) = - numnpnew
418  ENDIF
419  ENDDO
420  ENDDO
421  ptrtri3 => ptrtri3%next
422  ENDDO
423 
424  CALL com_new_attribute( surwin//'.nc', 'n', com_double, 3, 'm')
425  CALL com_set_size( surwin//'.nc', ibu, numnpnew )
426  CALL com_resize_array(surwin//'.nc', ibu, meshcoor, 3)
427 
428  DO i = 1, numnpnew
429  meshcoor(1:3,i) = coor(1:3,nodeflag2d(1,i))
430  END DO
431 
432 ! nodeflag to VolumeNode
433 
434  CALL com_new_attribute( surwin//'.bv', 'n', com_integer, 1, 'm')
435  CALL com_set_array( surwin//'.bv', ibu, nodeflag2d(2,1), 2)
436 
437 
438 ! connectiveity
439  testnumsurfel = numeltet2d(1,iprocs) ! 0
440 
441  CALL com_set_size( surwin//'.:t3', ibu, testnumsurfel)
442  CALL com_resize_array( surwin//'.:t3', ibu, elconntable, 3)
443 
444 
445 ! elemflag to VolumeElem
446 
447  CALL com_new_attribute( surwin//'.bf2c', 'e', com_integer, 1, '')
448  CALL com_resize_array(surwin//'.bf2c', ibu, elflag_list, 1)
449 
450  CALL com_new_attribute( surwin//'.bcflag', 'p', com_integer, 1, '')
451  CALL com_set_size( surwin//'.bcflag', ibu, 1)
452  CALL com_resize_array(surwin//'.bcflag', ibu, bcsurfflagone, 1)
453 
454  bcsurfflagone(1) = 1
455 
456  testnumsurfel = 0
457  ptrtri3 => surfmesh_tri3_sf_head
458  do while(associated(ptrtri3))
459  surfaceelonproc = .false.
460 !!$ do i = 1, 6
461 !!$ IF(NodeFlag(ptrtri6%ElemData(i)).eq.0) THEN
462 !!$ SurfaceElOnProc = .false.
463 !!$ exit
464 !!$ endif
465 !!$ enddo
466 
467  IF(epart(ptrtri3%ElemData(4)).EQ.iprocs) surfaceelonproc = .true.
468  IF(surfaceelonproc) THEN
469 !!$ WRITE(4001,'(10i10)') ABS(NodeFlag(ptrtri3%ElemData(1:3))),&
470 !!$ ElFlag(ptrtri3%ElemData(4)),1,1 ! last two are not really used
471 
472  testnumsurfel = testnumsurfel + 1
473  elconntable(1:3,testnumsurfel) = abs(nodeflag(ptrtri3%ElemData(1:3)))
474  elflag_list(testnumsurfel) = elflag(ptrtri3%ElemData(4))
475 
476 ! write(4003,*) ABS(NodeFlag(ptrtri3%ElemData(1))),ABS(NodeFlag(ptrtri3%ElemData(2))),ABS(NodeFlag(ptrtri3%ElemData(3)))
477  ENDIF
478  ptrtri3 => ptrtri3%next
479  ENDDO
480 
481  IF(testnumsurfel.NE.numeltet2d(1,iprocs))THEN
482  print*,'ERROR, number of triangles from link list in Ignitable surface mesh'
483  print*,' different then that in read_patran'
484  print*,' For Fluid Solid Inteface mesh'
485  print*,testnumsurfel,numeltet2d(1,iprocs)
486  stop
487  ENDIF
488 
489  nullify(ptrtri3)
490  ! obtain function handle ------------------------------------------------------
491 
492 ! write_attr = COM_get_function_handle( 'OUT.write_attribute')
493 ! set_option = COM_get_function_handle( 'OUT.set_option')
494 
495 ! IF(iProcs.EQ.1)THEN
496 
497  CALL com_call_function( set_option, 2, 'mode', 'w')
498 
499 ! ELSE
500 !
501 ! CALL COM_call_function( set_option, 2, 'mode', 'a')
502 !
503 ! ENDIF
504 
505 ! do not append process rank -----------------
506 
507 ! CALL COM_call_function( set_option, 2, 'rankwidth', '0')
508 
509 ! write surface window ------------------------
510 
511 ! sur_all = Com_get_attribute_handle( surWin//'.all')
512 ! CALL COM_call_function( write_attr, 4, 'Rocin/SurfMesh_FS_Ignt.hdf', sur_all,&
513 ! "solid_surf","00.000000")
514 
515 ! CALL COM_delete_window( surWin)
516 
517 
518 
519  ELSE IF(meshtype2d.EQ.4)THEN ! quad surface mesh ! 8 node brick
520 
521 
522 ! renumber elements nodes
523 
524  ptrhex8 => surfmesh_hex8_sf_head
525  DO WHILE(ASSOCIATED(ptrhex8).eqv..true.)
526 ! mark that changed surface mesh by changing it to a negative
527  DO i = 1, 4
528  glbndnum = ptrhex8%ElemData(i)
529  IF( nodeflag(glbndnum).GT.0) THEN ! if zero means not on this processor, neg means already renumber node
530  numnpnew = numnpnew + 1
531  nodeflag(glbndnum) = -nodeflag(glbndnum) ! overwrite nodeflag here but i need it for non S/F interface
532  ENDIF
533  ENDDO
534  ptrhex8 => ptrhex8%next
535  ENDDO
536 
537 ! write the surface mesh files
538 
539  IF(numelhex2d(1,iprocs).NE.0.AND.numnpnew.EQ.0)THEN
540  print*,'ERROR, found in surface mesh for Solid/Fluid Ignitable interface'
541  print*,'Found ',numelhex2d(1,iprocs),'quads but no nodes'
542  print*,'Stopping'
543  stop
544  ENDIF
545 
546  IF(numelhex2d(1,iprocs).EQ.0.AND.numnpnew.NE.0)THEN
547  print*,'Warning: found in sufrace mesh for Solid/Fluid Ignitable interface'
548  print*,'Found ',numnpnew,' nodes but no quads'
549  ENDIF
550 
551 
552 ! write the surface mesh files
553 
554  ! WRITE(4001,'(2i10)') NumNpNew, NumElhex2D(1,iProcs)
555 
556 !!$ IF(NumNPNew.NE.0.AND.NumElhex2D(1,iProcs).NE.0) &
557 !!$ write(4003,*) 'ZONE N=', NumNPNew, 'E=',NumElhex2D(1,iProcs),'F=FEPOINT, ET=QUADRILATERAL'
558 
559  nodeflag(:)= abs( nodeflag(:)) ! reset
560 
561  ALLOCATE(nodeflag2d(1:2,1:numnpnew)) ! for storing what I'm overwriting
562 
563  numnpnew = 0
564 ! renumber elements nodes
565 
566  testnumsurfel = 0
567 
568  ptrhex8 => surfmesh_hex8_sf_head
569  do while(associated(ptrhex8).eqv..true.)
570 ! mark that changed surface mesh by changing it to a negative
571  DO i = 1, 4
572  glbndnum = ptrhex8%ElemData(i)
573  IF( nodeflag(glbndnum).GT.0) THEN ! if zero means not on this processor
574  numnpnew = numnpnew + 1
575 ! WRITE(4001,'(3(x,e16.9),1i10)') coor(1:3,glbNdNum),NodeFlag(glbNdNum)
576 ! write(4003,*) coor(1:3,glbNdNum)
577  nodeflag2d(1,numnpnew) = glbndnum
578  nodeflag2d(2,numnpnew) = nodeflag(glbndnum)
579  nodeflag(glbndnum) = - numnpnew
580  endif
581  enddo
582  ptrhex8 => ptrhex8%next
583  ENDDO
584 
585  ! MS: Writing only non-vanishing windows
586  IF (numnpnew.GT.0) THEN
587  CALL com_new_attribute( surwin//'.nc', 'n', com_double, 3, 'm')
588  CALL com_set_size( surwin//'.nc', ibu, numnpnew )
589  CALL com_resize_array(surwin//'.nc', ibu, meshcoor, 3)
590 
591  DO i = 1, numnpnew
592  meshcoor(1:3,i) = coor(1:3,nodeflag2d(1,i))
593  END DO
594 
595  ! nodeflag to VolumeNode
596 
597  CALL com_new_attribute( surwin//'.bv', 'n', com_integer, 1, 'm')
598  CALL com_set_array( surwin//'.bv', ibu, nodeflag2d(2,1), 2)
599 
600  ! connectiveity
601 
602  testnumsurfel = numelhex2d(1,iprocs)
603 
604  CALL com_set_size( surwin//'.:q4', ibu, testnumsurfel)
605  CALL com_resize_array( surwin//'.:q4', ibu, elconntable, 4)
606 
607 
608  ! elemflag to VolumeElem
609 
610  CALL com_new_attribute( surwin//'.bf2c', 'e', com_integer, 1, '')
611  CALL com_resize_array(surwin//'.bf2c', ibu, elflag_list, 1)
612 
613  CALL com_new_attribute( surwin//'.bcflag', 'p', com_integer, 1, '')
614  CALL com_set_size( surwin//'.bcflag', ibu, 1)
615  CALL com_resize_array(surwin//'.bcflag', ibu, bcsurfflagone, 1)
616  bcsurfflagone(1) = 1
617  ENDIF
618  ! MS End
619 
620 
621  testnumsurfel = 0
622  ptrhex8 => surfmesh_hex8_sf_head
623  DO WHILE(ASSOCIATED(ptrhex8).eqv..true.)
624  surfaceelonproc = .false.
625 !!$ do i = 1, 4
626 !!$ IF(NodeFlag(ptrhex8%ElemData(i)).eq.0) THEN
627 !!$ SurfaceElOnProc = .false.
628 !!$ exit
629 !!$ endif
630 !!$ enddo
631  IF(epart(ptrhex8%ElemData(5)).EQ.iprocs) surfaceelonproc = .true.
632  IF(surfaceelonproc)THEN
633 ! WRITE(4001,'(10i10)') ABS(NodeFlag(ptrhex8%ElemData(1:4))),&
634 ! ElFlag(ptrhex8%ElemData(5)),1,1 ! last two are not really used
635 
636  testnumsurfel = testnumsurfel + 1
637  elconntable(1:4,testnumsurfel) = abs(nodeflag(ptrhex8%ElemData(1:4)))
638  elflag_list(testnumsurfel) = elflag(ptrhex8%ElemData(5))
639 ! write(4003,*) ABS(NodeFlag(ptrhex8%ElemData(1:4)))
640  ENDIF
641 
642  ptrhex8 => ptrhex8%next
643  ENDDO
644 
645  IF(testnumsurfel.NE.numelhex2d(1,iprocs))THEN
646  print*,'ERROR, number of quads from link list in Ignitable interface'
647  print*,' different then that in read_patran'
648  print*,testnumsurfel,numelhex2d(1,iprocs)
649  stop
650  ENDIF
651 
652  nullify(ptrhex8)
653 
654  CALL com_call_function( set_option, 2, 'mode', 'w')
655 
656  ENDIF
657  CALL com_call_function( set_option, 2, 'mode', 'w') ! why would this be called twice ? fix
658 
659  sur_all = com_get_attribute_handle( surwin//'.all')
660 
661 
662  ! MS: Writing only non-vanishing windows
663  IF (numnpnew.GT.0) THEN
664  IF (appendtoggle.EQV..true.) THEN
665  CALL com_call_function( set_option, 2, 'mode', 'a')
666  ELSE
667  CALL com_call_function( set_option, 2, 'mode', 'w')
668  appendtoggle = .true.
669  END IF
670  CALL com_call_function( write_attr, 4, 'Rocin/SurfMesh.'//ichr4, sur_all,&
671  "isolid","00.000000")
672  ENDIF
673  WRITE (*,*) 'iProc = ', iprocs, ' appendToggle = ', appendtoggle
674  ! MS : End
675 
676  CALL com_delete_window( surwin)
677 
678 ! --------------------------------------------------
679 ! ------------------------------------------------
680 ! (2) Non Ignitable interface (INTERACTING)
681 ! ---------------------------------------
682 ! -------------------------------------
683 
684  CALL com_new_window( surwin )
685 
686 
687 ! Reset NodeFlag for non solid/fluid interface
688 
689  IF(numnpnew.NE.0)THEN
690 
691  DO i = 1, numnpnew
692  nodeflag(nodeflag2d(1,i)) = nodeflag2d(2,i)
693  enddo
694  ENDIF
695  IF(associated(nodeflag2d).eqv..true.) deallocate(nodeflag2d)
696 
697 !!$ do i= 1, numnp_prmry
698 !!$ write(1001,*) i,NodeFlag(i)
699 !!$ enddo
700 
701  IF(meshtype2d.EQ.6)THEN
702 
703 ! renumber elements nodes
704 
705  ptrtri6 => surfmesh_tri6_sf_nonignt_head
706  DO WHILE(ASSOCIATED(ptrtri6))
707 ! mark that changed surface mesh by changing it to a negative
708  DO i = 1, 6
709  glbndnum = ptrtri6%ElemData(i)
710  IF( nodeflag(glbndnum).GT.0) THEN ! if zero means not on this processor, neg means already renumber node
711  numnpnew = numnpnew + 1
712  nodeflag(glbndnum) = -nodeflag(glbndnum) ! overwrite nodeflag here but i need it for non S/F interface
713  ENDIF
714  ENDDO
715  ptrtri6 => ptrtri6%next
716  ENDDO
717 
718 ! check for special situation where there is just an edge on the surface
719 
720  ptrtri6 => surfmesh_tri6_sf_nonignt_head
721  DO WHILE(ASSOCIATED(ptrtri6))
722  DO i = 1, 6
723  glbndnum = ptrtri6%ElemData(i)
724  DO k = 1,maxnumberofprocstosharenode
725  IF(procndlist(glbndnum,k).EQ.iprocs.AND.nodeflag(glbndnum).GT.0)THEN
726  print*,'Detected surface sliver for pressure'
727  numnpnew = numnpnew + 1
728  nodeflag(glbndnum) = -nodeflag(glbndnum) ! overwrite nodeflag here but i need it for non S/F interface
729  ENDIF
730  ENDDO
731  ENDDO
732  ptrtri6 => ptrtri6%next
733  ENDDO
734 
735 ! write the surface mesh files
736 
737  IF(numeltet2d(3,iprocs).NE.0.AND.numnpnew.EQ.0)THEN
738  print*,'ERROR, found in surface mesh for Solid/Fluid Ignitable interface'
739  print*,'Found ',numeltet2d(3,iprocs),'triangles but no nodes'
740  print*,'Stopping'
741  stop
742  ENDIF
743 
744  IF(numeltet2d(3,iprocs).EQ.0.AND.numnpnew.NE.0)THEN
745  print*,'Warning: found in sufrace mesh for Solid/Fluid Ignitable interface'
746  print*,'Found ',numnpnew,' nodes but no triangles'
747  ENDIF
748 
749 
750 
751 !!$ IF(NumNPNew.NE.0.AND.NumEltet2D(3,iProcs).NE.0) &
752 !!$ write(4003,*) 'ZONE N=', NumNPNew, 'E=',NumEltet2D(3,iProcs)*4,'F=FEPOINT, ET=TRIANGLE'
753 
754  nodeflag(:)= abs( nodeflag(:)) ! reset
755 
756  ALLOCATE(nodeflag2d(1:2,1:numnpnew)) ! for storing what I'm overwriting
757 
758  numnpnew = 0
759 ! renumber elements nodes
760 
761  ptrtri6 => surfmesh_tri6_sf_nonignt_head
762  DO WHILE(ASSOCIATED(ptrtri6).eqv..true.)
763 ! mark that changed surface mesh by changing it to a negative
764  DO i = 1, 6
765  glbndnum = ptrtri6%ElemData(i)
766  IF( nodeflag(glbndnum).GT.0) THEN ! if zero means not on this processor
767  numnpnew = numnpnew + 1
768  !WRITE(4001,'(3(1X,e16.9),1i10)') coor(1:3,glbNdNum),NodeFlag(glbNdNum)
769  !WRITE(4003,*) coor(1:3,glbNdNum)
770  nodeflag2d(1,numnpnew) = glbndnum
771  nodeflag2d(2,numnpnew) = nodeflag(glbndnum)
772  nodeflag(glbndnum) = - numnpnew
773  ENDIF
774  ENDDO
775  ptrtri6 => ptrtri6%next
776  ENDDO
777  ptrtri6 => surfmesh_tri6_sf_nonignt_head
778  DO WHILE(ASSOCIATED(ptrtri6))
779  DO i = 1, 6
780  glbndnum = ptrtri6%ElemData(i)
781  DO k = 1,maxnumberofprocstosharenode
782  IF(procndlist(glbndnum,k).EQ.iprocs.AND.nodeflag(glbndnum).GT.0)THEN
783  numnpnew = numnpnew + 1
784  print*,'Dectected surface sliver for pressure, writing'
785  !WRITE(4001,'(3(1X,e16.9),1i10)') coor(1:3,glbNdNum),NodeFlag(glbNdNum)
786  nodeflag2d(1,numnpnew) = glbndnum
787  nodeflag2d(2,numnpnew) = nodeflag(glbndnum)
788  nodeflag(glbndnum) = - numnpnew
789  ENDIF
790  ENDDO
791  ENDDO
792  ptrtri6 => ptrtri6%next
793  ENDDO
794 
795  CALL com_new_attribute( surwin//'.nc', 'n', com_double, 3, 'm')
796  CALL com_set_size( surwin//'.nc', inb, numnpnew )
797  CALL com_resize_array(surwin//'.nc', inb, meshcoor, 3)
798 
799  DO i = 1, numnpnew
800  meshcoor(1:3,i) = coor(1:3,nodeflag2d(1,i))
801  END DO
802 
803 
804  CALL com_new_attribute( surwin//'.bv', 'n', com_integer, 1, '')
805 ! CALL COM_set_size( surWin//'.bv', iNB, NumNpNew )
806  CALL com_set_array(surwin//'.bv', inb, nodeflag2d(2,1), 2)
807 
808 ! connectiveity
809 
810 
811  testnumsurfel = numeltet2d(3,iprocs)
812 
813  CALL com_set_size( surwin//'.:t6', inb, testnumsurfel)
814  CALL com_resize_array( surwin//'.:t6', inb, elconntable, 6)
815 
816 
817  CALL com_new_attribute( surwin//'.bf2c', 'e', com_integer, 1, '')
818  CALL com_resize_array(surwin//'.bf2c', inb, elflag_list, 1)
819 
820  CALL com_new_attribute( surwin//'.bcflag', 'p', com_integer, 1, '')
821  CALL com_set_size( surwin//'.bcflag', inb, 1)
822  CALL com_resize_array(surwin//'.bcflag', inb, bcsurfflagzero, 1)
823 
824  bcsurfflagzero(1) = 0
825 
826 
827 
828  testnumsurfel = 0
829  ptrtri6 => surfmesh_tri6_sf_nonignt_head
830  DO WHILE(ASSOCIATED(ptrtri6))
831  surfaceelonproc = .false.
832 
833  IF(epart(ptrtri6%ElemData(7)).EQ.iprocs) surfaceelonproc = .true.
834  IF(surfaceelonproc) THEN
835 !!$ WRITE(4001,'(10i10)') ABS(NodeFlag(ptrtri6%ElemData(1:6))),&
836 !!$ ElFlag(ptrtri6%ElemData(7)),1,1 ! last two are not really used
837 
838  testnumsurfel = testnumsurfel + 1
839  elconntable(1:6,testnumsurfel) = abs(nodeflag(ptrtri6%ElemData(1:6)))
840  elflag_list(testnumsurfel) = elflag(ptrtri6%ElemData(7))
841  !write(4003,*) ABS(NodeFlag(ptrtri6%ElemData(1))),ABS(NodeFlag(ptrtri6%ElemData(4))),ABS(NodeFlag(ptrtri6%ElemData(6)))
842  !write(4003,*) ABS(NodeFlag(ptrtri6%ElemData(4))),ABS(NodeFlag(ptrtri6%ElemData(2))),ABS(NodeFlag(ptrtri6%ElemData(5)))
843  !write(4003,*) ABS(NodeFlag(ptrtri6%ElemData(6))),ABS(NodeFlag(ptrtri6%ElemData(4))),ABS(NodeFlag(ptrtri6%ElemData(5)))
844  !write(4003,*) ABS(NodeFlag(ptrtri6%ElemData(3))),ABS(NodeFlag(ptrtri6%ElemData(6))),ABS(NodeFlag(ptrtri6%ElemData(5)))
845  ENDIF
846  ptrtri6 => ptrtri6%next
847  ENDDO
848 
849  IF(testnumsurfel.NE.numeltet2d(3,iprocs))THEN
850  print*,'ERROR, number of triangles from link list in Mesh2d'
851  print*,' different then that in read_patran'
852  print*,testnumsurfel,numeltet2d(3,iprocs)
853  stop
854  ENDIF
855 
856  nullify(ptrtri6)
857 
858 
859  ! obtain function handle ------------------------------------------------------
860 
861 ! write_attr = COM_get_function_handle( 'OUT.write_attribute')
862 ! set_option = COM_get_function_handle( 'OUT.set_option')
863 
864 ! IF(iProcs.EQ.1)THEN
865 
866 ! CALL COM_call_function( set_option, 2, 'mode', 'w')
867 
868 ! ELSE
869 
870 ! CALL COM_call_function( set_option, 2, 'mode', 'a')
871 
872 ! ENDIF
873 
874 ! do not append process rank -----------------
875 
876 ! CALL COM_call_function( set_option, 2, 'rankwidth', '0')
877 
878 ! write surface window ------------------------
879 
880 ! sur_all = Com_get_attribute_handle( surWin//'.all')
881 
882 ! CALL COM_call_function( write_attr, 4, 'Rocin/SurfMesh_FS_NonIgnt.hdf', sur_all,&
883 ! "solid_surf","00.000000")
884 
885 ! CALL COM_delete_window( surWin)
886 
887 
888  ELSE IF(meshtype2d.EQ.3)THEN
889 
890 ! renumber elements nodes
891 
892  ptrtri3 => surfmesh_tri3_sf_nonignt_head
893  do while(associated(ptrtri3).eqv..true.)
894 ! mark that changed surface mesh by changing it to a negative
895  DO i = 1, 3
896  glbndnum = ptrtri3%ElemData(i)
897  IF( nodeflag(glbndnum).GT.0) THEN ! if zero means not on this processor, neg means already renumber node
898  numnpnew = numnpnew + 1
899  nodeflag(glbndnum) = -nodeflag(glbndnum) ! overwrite nodeflag here but i need it for non S/F interface
900  endif
901  enddo
902  ptrtri3 => ptrtri3%next
903  ENDDO
904 
905 ! check for special situation where there is just an edge on the surface
906 
907  ptrtri3 => surfmesh_tri3_sf_nonignt_head
908  do while(associated(ptrtri3).eqv..true.)
909  DO i = 1, 3
910  glbndnum = ptrtri3%ElemData(i)
911  DO k = 1,maxnumberofprocstosharenode
912  IF(procndlist(glbndnum,k).EQ.iprocs.AND.nodeflag(glbndnum).GT.0)THEN
913  print*,'Dectected surface sliver for pressure'
914  numnpnew = numnpnew + 1
915  nodeflag(glbndnum) = -nodeflag(glbndnum) ! overwrite nodeflag here but i need it for non S/F interface
916 
917  ENDIF
918  ENDDO
919  ENDDO
920  ptrtri3 => ptrtri3%next
921  ENDDO
922 
923 ! write the surface mesh files
924 
925 
926  !WRITE(4001,'(2i10)') NumNpNew, NumEltet2D(3,iProcs)
927 
928  IF(numeltet2d(3,iprocs).NE.0.AND.numnpnew.EQ.0)THEN
929  print*,'ERROR, found in surface mesh for Solid/Fluid interface'
930  print*,'Found ',numeltet2d(3,iprocs),'triangles but no nodes'
931  print*,'Stopping'
932  stop
933  ENDIF
934 
935  IF(numeltet2d(3,iprocs).EQ.0.AND.numnpnew.NE.0)THEN
936  print*,'Warning: found in sufrace mesh for Solid/Fluid interface'
937  print*,'Found ',numnpnew,' nodes but no triangles'
938  ENDIF
939 
940 ! IF(NumNPNew.NE.0.AND.NumEltet2D(3,iProcs).NE.0) &
941 ! write(4003,*) 'ZONE N=', NumNPNew, 'E=',NumEltet2D(3,iProcs),'F=FEPOINT, ET=TRIANGLE'
942 
943  nodeflag(:)= abs( nodeflag(:)) ! reset
944 
945  ALLOCATE(nodeflag2d(1:2,1:numnpnew)) ! for storing what I'm overwriting
946 
947  numnpnew = 0
948 ! renumber elements nodes
949 
950  ptrtri3 => surfmesh_tri3_sf_nonignt_head
951  do while(associated(ptrtri3).eqv..true.)
952 ! mark that changed surface mesh by changing it to a negative
953  DO i = 1, 3
954  glbndnum = ptrtri3%ElemData(i)
955  IF( nodeflag(glbndnum).GT.0) THEN ! if zero means not on this processor
956  numnpnew = numnpnew + 1
957 ! WRITE(4001,'(3(x,e16.9),1i10)') coor(1:3,glbNdNum),NodeFlag(glbNdNum)
958 ! write(4003,*) coor(1:3,glbNdNum)
959  nodeflag2d(1,numnpnew) = glbndnum
960  nodeflag2d(2,numnpnew) = nodeflag(glbndnum)
961  nodeflag(glbndnum) = - numnpnew
962  endif
963  enddo
964  ptrtri3 => ptrtri3%next
965  ENDDO
966  ptrtri3 => surfmesh_tri3_sf_nonignt_head
967  do while(associated(ptrtri3).eqv..true.)
968  DO i = 1, 3
969  glbndnum = ptrtri3%ElemData(i)
970  DO k = 1,maxnumberofprocstosharenode
971  IF(procndlist(glbndnum,k).EQ.iprocs.AND.nodeflag(glbndnum).GT.0)THEN
972  numnpnew = numnpnew + 1
973  print*,'Dectected surface sliver for pressure, writing'
974 ! WRITE(4001,'(3(x,e16.9),1i10)') coor(1:3,glbNdNum),NodeFlag(glbNdNum)
975  nodeflag2d(1,numnpnew) = glbndnum
976  nodeflag2d(2,numnpnew) = nodeflag(glbndnum)
977  nodeflag(glbndnum) = - numnpnew
978  ENDIF
979  ENDDO
980  ENDDO
981  ptrtri3 => ptrtri3%next
982  ENDDO
983 
984  CALL com_new_attribute( surwin//'.nc', 'n', com_double, 3, 'm')
985  CALL com_set_size( surwin//'.nc', inb, numnpnew )
986  CALL com_resize_array(surwin//'.nc', inb, meshcoor, 3)
987 
988  DO i = 1, numnpnew
989  meshcoor(1:3,i) = coor(1:3,nodeflag2d(1,i))
990  END DO
991 
992 
993  CALL com_new_attribute( surwin//'.bv', 'n', com_integer, 1, '')
994 ! CALL COM_set_size( surWin//'.bv', iNB, NumNpNew )
995  CALL com_set_array(surwin//'.bv', inb, nodeflag2d(2,1), 2)
996 ! connectiveity
997 
998  testnumsurfel = numeltet2d(3,iprocs)
999 
1000  CALL com_set_size( surwin//'.:t3', inb, testnumsurfel)
1001  CALL com_resize_array( surwin//'.:t3', inb, elconntable, 3)
1002 
1003 
1004  CALL com_new_attribute( surwin//'.bf2c', 'e', com_integer, 1, '')
1005  CALL com_resize_array(surwin//'.bf2c', inb, elflag_list, 1)
1006 
1007  CALL com_new_attribute( surwin//'.bcflag', 'p', com_integer, 1, '')
1008  CALL com_set_size( surwin//'.bcflag', inb, 1)
1009  CALL com_resize_array(surwin//'.bcflag', inb, bcsurfflagzero, 1)
1010 
1011  bcsurfflagzero(1) = 0
1012 
1013  testnumsurfel = 0
1014  ptrtri3 => surfmesh_tri3_sf_nonignt_head
1015  do while(associated(ptrtri3).eqv..true.)
1016  surfaceelonproc = .false.
1017 !!$ do i = 1, 6
1018 !!$ IF(NodeFlag(ptrtri6%ElemData(i)).eq.0) THEN
1019 !!$ SurfaceElOnProc = .false.
1020 !!$ exit
1021 !!$ endif
1022 !!$ enddo
1023 
1024  IF(epart(ptrtri3%ElemData(4)).EQ.iprocs) surfaceelonproc = .true.
1025  IF(surfaceelonproc.eqv..true.) THEN
1026 ! WRITE(4001,'(10i10)') ABS(NodeFlag(ptrtri3%ElemData(1:3))),&
1027 ! ElFlag(ptrtri3%ElemData(4)),1,1 ! last two are not really used
1028  testnumsurfel = testnumsurfel + 1
1029  elconntable(1:3,testnumsurfel) = abs(nodeflag(ptrtri3%ElemData(1:3)))
1030  elflag_list(testnumsurfel) = elflag(ptrtri3%ElemData(4))
1031  ENDIF
1032  ptrtri3 => ptrtri3%next
1033  ENDDO
1034 
1035  IF(testnumsurfel.NE.numeltet2d(3,iprocs))THEN
1036  print*,'ERROR, number of triangles from link list in Mesh2d'
1037  print*,' different then that in read_patran'
1038  print*,' For Fluid Solid Inteface mesh'
1039  print*,testnumsurfel,numeltet2d(3,iprocs)
1040  stop
1041  ENDIF
1042 
1043  nullify(ptrtri3)
1044 
1045  ! obtain function handle ------------------------------------------------------
1046 
1047 ! write_attr = COM_get_function_handle( 'OUT.write_attribute')
1048 ! set_option = COM_get_function_handle( 'OUT.set_option')
1049 
1050 ! IF(iProcs.EQ.1)THEN
1051 
1052 ! CALL COM_call_function( set_option, 2, 'mode', 'w')
1053 
1054 ! ELSE
1055 
1056 ! CALL COM_call_function( set_option, 2, 'mode', 'a')
1057 
1058 ! ENDIF
1059 
1060 ! do not append process rank -----------------
1061 
1062 ! CALL COM_call_function( set_option, 2, 'rankwidth', '0')
1063 
1064 ! write surface window ------------------------
1065 
1066 ! sur_all = Com_get_attribute_handle( surWin//'.all')
1067 
1068 ! CALL COM_call_function( write_attr, 4, 'Rocin/SurfMesh_FS_NonIgnt.hdf', sur_all,&
1069 ! "solid_surf","00.000000")
1070 
1071 ! CALL COM_delete_window( surWin)
1072 
1073  ELSE IF(meshtype2d.EQ.4)THEN ! quad surface mesh ! 8 node brick
1074 
1075 
1076 ! renumber elements nodes
1077 
1078  ptrhex8 => surfmesh_hex8_sf_nonignt_head
1079  DO WHILE(ASSOCIATED(ptrhex8).eqv..true.)
1080 ! mark that changed surface mesh by changing it to a negative
1081  DO i = 1, 4
1082  glbndnum = ptrhex8%ElemData(i)
1083  IF( nodeflag(glbndnum).GT.0) THEN ! if zero means not on this processor, neg means already renumber node
1084  numnpnew = numnpnew + 1
1085  nodeflag(glbndnum) = -nodeflag(glbndnum) ! overwrite nodeflag here but i need it for non S/F interface
1086  ENDIF
1087  ENDDO
1088  ptrhex8 => ptrhex8%next
1089  ENDDO
1090 
1091 ! write the surface mesh files
1092 
1093  IF(numelhex2d(3,iprocs).NE.0.AND.numnpnew.EQ.0)THEN
1094  print*,'ERROR, found in surface mesh for Solid/Fluid Ignitable interface'
1095  print*,'Found ',numelhex2d(3,iprocs),'quads but no nodes'
1096  print*,'Stopping'
1097  stop
1098  ENDIF
1099 
1100  IF(numelhex2d(3,iprocs).EQ.0.AND.numnpnew.NE.0)THEN
1101  print*,'Warning: found in sufrace mesh for Solid/Fluid Ignitable interface'
1102  print*,'Found ',numnpnew,' nodes but no quads'
1103  ENDIF
1104 
1105 ! WRITE(4001,'(2i10)') NumNpNew, NumElhex2D(3,iProcs)
1106 
1107 ! IF(NumNPNew.NE.0.AND.NumElhex2D(3,iProcs).NE.0) &
1108 ! write(4003,*) 'ZONE N=', NumNPNew, 'E=',NumElhex2D(3,iProcs),'F=FEPOINT, ET=QUADRILATERAL'
1109 
1110  nodeflag(:)= abs( nodeflag(:)) ! reset
1111 
1112  ALLOCATE(nodeflag2d(1:2,1:numnpnew)) ! for storing what I'm overwriting
1113 
1114  numnpnew = 0
1115 ! renumber elements nodes
1116  ptrhex8 => surfmesh_hex8_sf_nonignt_head
1117  do while(associated(ptrhex8).eqv..true.)
1118 ! mark that changed surface mesh by changing it to a negative
1119  DO i = 1, 4
1120  glbndnum = ptrhex8%ElemData(i)
1121  IF( nodeflag(glbndnum).GT.0) THEN ! if zero means not on this processor
1122  numnpnew = numnpnew + 1
1123 ! WRITE(4001,'(3(x,e16.9),1i10)') coor(1:3,glbNdNum),NodeFlag(glbNdNum)
1124 ! write(4003,*) coor(1:3,glbNdNum)
1125  nodeflag2d(1,numnpnew) = glbndnum
1126  nodeflag2d(2,numnpnew) = nodeflag(glbndnum)
1127  nodeflag(glbndnum) = - numnpnew
1128  endif
1129  enddo
1130  ptrhex8 => ptrhex8%next
1131  ENDDO
1132 
1133  ! MS
1134  IF (numnpnew.GT.0) THEN
1135  CALL com_new_attribute( surwin//'.nc', 'n', com_double, 3, 'm')
1136  CALL com_set_size( surwin//'.nc', inb, numnpnew )
1137  CALL com_resize_array(surwin//'.nc', inb, meshcoor, 3)
1138 
1139  DO i = 1, numnpnew
1140  meshcoor(1:3,i) = coor(1:3,nodeflag2d(1,i))
1141  END DO
1142 
1143  CALL com_new_attribute( surwin//'.bv', 'n', com_integer, 1, '')
1144  ! CALL COM_set_size( surWin//'.bv', iNB, NumNpNew )
1145  CALL com_set_array(surwin//'.bv', inb, nodeflag2d(2,1), 2)
1146 
1147  ! connectiveity
1148 
1149 
1150  testnumsurfel = numelhex2d(3,iprocs)
1151 
1152  CALL com_set_size( surwin//'.:q4', inb, testnumsurfel)
1153  CALL com_resize_array( surwin//'.:q4', inb, elconntable, 4)
1154 
1155 
1156  CALL com_new_attribute( surwin//'.bf2c', 'e', com_integer, 1, '')
1157  CALL com_resize_array(surwin//'.bf2c', inb, elflag_list, 1)
1158 
1159  CALL com_new_attribute( surwin//'.bcflag', 'p', com_integer, 1, '')
1160  CALL com_set_size( surwin//'.bcflag', inb, 1)
1161  CALL com_resize_array(surwin//'.bcflag', inb, bcsurfflagzero, 1)
1162  bcsurfflagzero(1) = 0
1163  END IF
1164  ! MS End
1165 
1166  testnumsurfel = 0
1167  ptrhex8 => surfmesh_hex8_sf_nonignt_head
1168  DO WHILE(ASSOCIATED(ptrhex8).eqv..true.)
1169  surfaceelonproc = .false.
1170 
1171  IF(epart(ptrhex8%ElemData(5)).EQ.iprocs) surfaceelonproc = .true.
1172  IF(surfaceelonproc.eqv..true.)THEN
1173 ! WRITE(4001,'(10i10)') ABS(NodeFlag(ptrhex8%ElemData(1:4))),&
1174 ! ElFlag(ptrhex8%ElemData(5)),1,1 ! last two are not really used
1175 
1176  testnumsurfel = testnumsurfel + 1
1177  elconntable(1:4,testnumsurfel) = abs(nodeflag(ptrhex8%ElemData(1:4)))
1178  elflag_list(testnumsurfel) = elflag(ptrhex8%ElemData(5))
1179 
1180 ! write(4003,*) ABS(NodeFlag(ptrhex8%ElemData(1:4)))
1181  ENDIF
1182  ptrhex8 => ptrhex8%next
1183  ENDDO
1184 
1185  IF(testnumsurfel.NE.numelhex2d(3,iprocs))THEN
1186  print*,'ERROR, number of quads from link list in Non-Ignitable surface mesh'
1187  print*,' different then that in read_patran'
1188  print*,testnumsurfel,numelhex2d(3,iprocs)
1189  stop
1190  ENDIF
1191 
1192  nullify(ptrhex8)
1193 
1194  ENDIF
1195 
1196 ! Reset NodeFlag for non solid/fluid interface
1197 
1198  IF(numnpnew.NE.0)THEN
1199 
1200  DO i = 1, numnpnew
1201  nodeflag(nodeflag2d(1,i)) = nodeflag2d(2,i)
1202  enddo
1203  ENDIF
1204 
1205  ! write surface window ------------------------
1206 
1207  CALL com_call_function( set_option, 2, 'mode', 'a')
1208 
1209  sur_all = com_get_attribute_handle( surwin//'.all')
1210 
1211  ! MS
1212  IF (numnpnew.GT.0) THEN
1213  IF (appendtoggle.EQV..true.) THEN
1214  CALL com_call_function( set_option, 2, 'mode', 'a')
1215  ELSE
1216  CALL com_call_function( set_option, 2, 'mode', 'w')
1217  appendtoggle = .true.
1218  END IF
1219  CALL com_call_function( write_attr, 4, 'Rocin/SurfMesh.'//ichr4, sur_all,&
1220  "isolid","00.000000")
1221  END IF
1222  WRITE (*,*) 'iProc = ', iprocs, ' appendToggle = ', appendtoggle
1223  ! MS End
1224 
1225  CALL com_delete_window( surwin)
1226 
1227  IF(associated(nodeflag2d).eqv..true.) deallocate(nodeflag2d)
1228 
1229 ! --------------------------------------------------
1230 ! ------------------------------------------------
1231 ! (3) Non Fluid/Structure interface
1232 ! ---------------------------------------
1233 ! -------------------------------------
1234 
1235 ! renumber elements nodes
1236 
1237  CALL com_new_window( surwin)
1238 
1239  IF(meshtype2d.EQ.6)THEN
1240 
1241  numnpnew = 0
1242  nullify(ptrtri6)
1243 ! renumber elements nodes
1244  ptrtri6 => surfmesh_tri6_s_head
1245 
1246  DO WHILE(ASSOCIATED(ptrtri6).eqv..true.)
1247 ! mark that changed surface mesh by changing it to a negative
1248  DO i = 1, 6
1249  glbndnum = ptrtri6%ElemData(i)
1250  IF( nodeflag(glbndnum).GT.0) THEN ! if zero means not on this processor
1251  numnpnew = numnpnew + 1
1252  nodeflag(glbndnum) = -nodeflag(glbndnum)
1253  ENDIF
1254  ENDDO
1255  ptrtri6 => ptrtri6%next
1256  ENDDO
1257 
1258 
1259 ! check for special situation where there is just an edge on the surface
1260 
1261  ptrtri6 => surfmesh_tri6_s_head
1262  do while(associated(ptrtri6).eqv..true.)
1263  DO i = 1, 6
1264  glbndnum = ptrtri6%ElemData(i)
1265  DO k = 1,maxnumberofprocstosharenode
1266  IF(procndlist(glbndnum,k).EQ.iprocs.AND.nodeflag(glbndnum).GT.0)THEN
1267  print*,'Dectected surface sliver for pressure'
1268  numnpnew = numnpnew + 1
1269  nodeflag(glbndnum) = -nodeflag(glbndnum)
1270  ENDIF
1271  ENDDO
1272  ENDDO
1273  ptrtri6 => ptrtri6%next
1274  ENDDO
1275 
1276 ! WRITE(4002,'(2i10)') NumNpNew, NumEltet2D(2,iProcs)
1277 ! write(4003,*) 'ZONE N=', NumNPNew, 'E=',NumEltet2D(2,iProcs),'F=FEPOINT, ET=TRIANGLE'
1278 
1279  nodeflag(:)= abs( nodeflag(:)) ! reset
1280 
1281  ALLOCATE(nodeflag2d(1:2,1:numnpnew)) ! for storing what I'm overwriting
1282 
1283  numnpnew = 0
1284 ! renumber elements nodes
1285 
1286  ptrtri6 => surfmesh_tri6_s_head
1287 
1288  do while(associated(ptrtri6).eqv..true.)
1289 ! mark that changed surface mesh by changing it to a negative
1290  DO i = 1, 6
1291  glbndnum = ptrtri6%ElemData(i)
1292  IF( nodeflag(glbndnum).GT.0) THEN ! if zero means not on this processor
1293  numnpnew = numnpnew + 1
1294 ! WRITE(4002,'(3(x,e16.9),1i10)') coor(1:3,glbNdNum),NodeFlag(glbNdNum)
1295 ! write(4003,*) coor(1:3,glbNdNum)
1296  nodeflag2d(1,numnpnew) = glbndnum
1297  nodeflag2d(2,numnpnew) = nodeflag(glbndnum)
1298  nodeflag(ptrtri6%ElemData(i)) = - numnpnew
1299  endif
1300  enddo
1301  ptrtri6 => ptrtri6%next
1302  ENDDO
1303  ptrtri6 => surfmesh_tri6_s_head
1304  do while(associated(ptrtri6).eqv..true.)
1305  DO i = 1, 6
1306  glbndnum = ptrtri6%ElemData(i)
1307  DO k = 1,maxnumberofprocstosharenode
1308  IF(procndlist(glbndnum,k).EQ.iprocs.AND.nodeflag(glbndnum).GT.0)THEN
1309  print*,'Dectected surface sliver for pressure'
1310  numnpnew = numnpnew + 1
1311 ! WRITE(4002,'(3(x,e16.9),1i10)') coor(1:3,glbNdNum),NodeFlag(glbNdNum)
1312  nodeflag2d(1,numnpnew) = glbndnum
1313  nodeflag2d(2,numnpnew) = nodeflag(glbndnum)
1314  nodeflag(glbndnum) = - numnpnew
1315  ENDIF
1316  ENDDO
1317  ENDDO
1318  ptrtri6 => ptrtri6%next
1319  ENDDO
1320 
1321 
1322  CALL com_new_attribute( surwin//'.nc', 'n', com_double, 3, 'm')
1323  CALL com_set_size( surwin//'.nc', ini, numnpnew )
1324  CALL com_resize_array(surwin//'.nc', ini, meshcoor, 3)
1325 
1326  DO i = 1, numnpnew
1327  meshcoor(1:3,i) = coor(1:3,nodeflag2d(1,i))
1328  END DO
1329 
1330 
1331  CALL com_new_attribute( surwin//'.bv', 'n', com_integer, 1, '')
1332 ! CALL COM_set_size( surWin//'.bv', iNI, NumNpNew )
1333  CALL com_set_array( surwin//'.bv', ini, nodeflag2d(2,1), 2)
1334 
1335 
1336 ! connectiveity
1337  testnumsurfel = numeltet2d(2,iprocs)
1338  print*,'Number of Non-Interacting nodes and elements',numnpnew,testnumsurfel
1339 
1340  CALL com_set_size( surwin//'.:t6', ini, testnumsurfel)
1341  CALL com_resize_array( surwin//'.:t6', ini, elconntable, 6)
1342 
1343 
1344  CALL com_new_attribute( surwin//'.bf2c', 'e', com_integer, 1, '')
1345  CALL com_resize_array(surwin//'.bf2c', ini, elflag_list, 1)
1346 
1347  CALL com_new_attribute( surwin//'.bcflag', 'p', com_integer, 1, '')
1348  CALL com_set_size( surwin//'.bcflag', ini, 1)
1349  CALL com_resize_array(surwin//'.bcflag', ini, bcsurfflagtwo, 1)
1350 
1351  bcsurfflagtwo(1) = 2
1352 
1353  testnumsurfel = 0
1354  ptrtri6 => surfmesh_tri6_s_head
1355  do while(associated(ptrtri6).eqv..true.)
1356  surfaceelonproc = .false.
1357 !!$ do i = 1, 6
1358 !!$ IF(NodeFlag(ptrtri6%ElemData(i)).eq.0)THEN
1359 !!$ SurfaceElOnProc = .false.
1360 !!$ exit
1361 !!$ endif
1362 !!$ enddo
1363  IF(epart(ptrtri6%ElemData(7)).EQ.iprocs) surfaceelonproc = .true.
1364  IF(surfaceelonproc) THEN
1365  ! WRITE(4002,'(10i10)') ABS(NodeFlag(ptrtri6%ElemData(1:6))),&
1366 ! ElFlag(ptrtri6%ElemData(7)),2,1 ! last two are not really used
1367 !!$ IF(SurfaceElOnProc) THEN
1368 !!$ write(4003,*) ABS(NodeFlag(ptrtri6%ElemData(1:6)))
1369 !!$ ENDIF
1370  testnumsurfel = testnumsurfel + 1
1371  elconntable(1:6,testnumsurfel) = abs(nodeflag(ptrtri6%ElemData(1:6)))
1372  ENDIF
1373  ptrtri6 => ptrtri6%next
1374  ENDDO
1375 
1376  IF(testnumsurfel.NE.numeltet2d(2,iprocs))THEN
1377  print*,'ERROR, number of triangles from link list in Mesh2d'
1378  print*,' different then that in read_patran'
1379  print*,testnumsurfel,numeltet2d(2,iprocs)
1380  stop
1381  ENDIF
1382 
1383  nullify(ptrtri6) ! obtain function handle ------------------------------------------------------
1384 
1385 ! write_attr = COM_get_function_handle( 'OUT.write_attribute')
1386 ! set_option = COM_get_function_handle( 'OUT.set_option')
1387 
1388 ! IF(iProcs.EQ.1)THEN
1389 
1390 ! CALL COM_call_function( set_option, 2, 'mode', 'w')
1391 
1392 ! ELSE
1393 
1394 ! CALL COM_call_function( set_option, 2, 'mode', 'a')
1395 
1396 ! ENDIF
1397 
1398 ! do not append process rank -----------------
1399 
1400 ! CALL COM_call_function( set_option, 2, 'rankwidth', '0')
1401 
1402 ! write surface window ------------------------
1403 
1404 ! sur_all = Com_get_attribute_handle( surWin//'.all')
1405 ! CALL COM_call_function( write_attr, 4, 'Rocin/SurfMesh_S.hdf', sur_all,&
1406 ! "solid_surf","00.000000")
1407 
1408 
1409 
1410  ELSE IF(meshtype2d.EQ.3)THEN
1411 
1412  numnpnew = 0
1413  nullify(ptrtri3)
1414 ! renumber elements nodes
1415 
1416  ptrtri3 => surfmesh_tri3_s_head
1417 
1418  do while(associated(ptrtri3).eqv..true.)
1419 ! mark that changed surface mesh by changing it to a negative
1420  DO i = 1, 3
1421  glbndnum = ptrtri3%ElemData(i)
1422  IF( nodeflag(glbndnum).GT.0) THEN ! if zero means not on this processor
1423  numnpnew = numnpnew + 1
1424  nodeflag(glbndnum) = -nodeflag(glbndnum)
1425  endif
1426  enddo
1427  ptrtri3 => ptrtri3%next
1428  ENDDO
1429 
1430 
1431 ! check for special situation where there is just an edge on the surface
1432 
1433  ptrtri3 => surfmesh_tri3_s_head
1434  do while(associated(ptrtri3).eqv..true.)
1435  DO i = 1, 3
1436  glbndnum = ptrtri3%ElemData(i)
1437  DO k = 1,maxnumberofprocstosharenode
1438  IF(procndlist(glbndnum,k).EQ.iprocs.AND.nodeflag(glbndnum).GT.0)THEN
1439  print*,'Dectected surface sliver for pressure'
1440  numnpnew = numnpnew + 1
1441  nodeflag(glbndnum) = -nodeflag(glbndnum)
1442  ENDIF
1443  ENDDO
1444  ENDDO
1445  ptrtri3 => ptrtri3%next
1446  ENDDO
1447 ! WRITE(4002,'(2i10)') NumNpNew, NumEltet2D(2,iProcs)
1448 ! write(4003,*) 'ZONE N=', NumNPNew, 'E=',NumEltet2D(2,iProcs),'F=FEPOINT, ET=TRIANGLE'
1449 
1450  nodeflag(:)= abs( nodeflag(:)) ! reset
1451 
1452  ALLOCATE(nodeflag2d(1:2,1:numnpnew)) ! for storing what I'm overwriting
1453 
1454  numnpnew = 0
1455 ! renumber elements nodes
1456 
1457  ptrtri3 => surfmesh_tri3_s_head
1458 
1459  do while(associated(ptrtri3).eqv..true.)
1460 ! mark that changed surface mesh by changing it to a negative
1461  DO i = 1, 3
1462  glbndnum = ptrtri3%ElemData(i)
1463  IF( nodeflag(glbndnum).GT.0) THEN ! if zero means not on this processor
1464  numnpnew = numnpnew + 1
1465 ! WRITE(4002,'(3(x,e16.9),1i10)') coor(1:3,glbNdNum),NodeFlag(glbNdNum)
1466 ! write(4003,*) coor(1:3,glbNdNum)
1467  nodeflag2d(1,numnpnew) = glbndnum
1468  nodeflag2d(2,numnpnew) = nodeflag(glbndnum)
1469  nodeflag(ptrtri3%ElemData(i)) = - numnpnew
1470  endif
1471  enddo
1472  ptrtri3 => ptrtri3%next
1473  ENDDO
1474  ptrtri3 => surfmesh_tri3_s_head
1475  do while(associated(ptrtri3).eqv..true.)
1476  DO i = 1, 3
1477  glbndnum = ptrtri3%ElemData(i)
1478  DO k = 1,maxnumberofprocstosharenode
1479  IF(procndlist(glbndnum,k).EQ.iprocs.AND.nodeflag(glbndnum).GT.0)THEN
1480  print*,'Dectected surface sliver for pressure'
1481  numnpnew = numnpnew + 1
1482 ! WRITE(4002,'(3(x,e16.9),1i10)') coor(1:3,glbNdNum),NodeFlag(glbNdNum)
1483  nodeflag2d(1,numnpnew) = glbndnum
1484  nodeflag2d(2,numnpnew) = nodeflag(glbndnum)
1485  nodeflag(ptrtri3%ElemData(i)) = - numnpnew
1486  ENDIF
1487  ENDDO
1488  ENDDO
1489  ptrtri3 => ptrtri3%next
1490  ENDDO
1491 
1492 
1493  CALL com_new_attribute( surwin//'.nc', 'n', com_double, 3, 'm')
1494  CALL com_set_size( surwin//'.nc', ini, numnpnew )
1495  CALL com_resize_array(surwin//'.nc', ini, meshcoor, 3)
1496 
1497  DO i = 1, numnpnew
1498  meshcoor(1:3,i) = coor(1:3,nodeflag2d(1,i))
1499  END DO
1500 
1501  CALL com_new_attribute( surwin//'.bv', 'n', com_integer, 1, '')
1502 ! CALL COM_set_size( surWin//'.bv', iNI, NumNpNew )
1503  CALL com_set_array( surwin//'.bv', ini, nodeflag2d(2,1), 2)
1504 
1505 
1506 ! connectiveity
1507  testnumsurfel = numeltet2d(2,iprocs) !0
1508 
1509  CALL com_set_size( surwin//'.:t3', ini, testnumsurfel)
1510  CALL com_resize_array( surwin//'.:t3', ini, elconntable, 3)
1511 
1512 
1513  CALL com_new_attribute( surwin//'.bf2c', 'e', com_integer, 1, '')
1514  CALL com_resize_array(surwin//'.bf2c', ini, elflag_list, 1)
1515 
1516 
1517  CALL com_new_attribute( surwin//'.bcflag', 'p', com_integer, 1, '')
1518  CALL com_set_size( surwin//'.bcflag', ini,1)
1519  CALL com_resize_array(surwin//'.bcflag', ini, bcsurfflagtwo, 1)
1520 
1521  bcsurfflagtwo(1) = 2
1522 
1523  testnumsurfel = 0
1524  ptrtri3 => surfmesh_tri3_s_head
1525  do while(associated(ptrtri3).eqv..true.)
1526  surfaceelonproc = .false.
1527 !!$ do i = 1, 3
1528 !!$ IF(NodeFlag(ptrtri3%ElemData(i)).eq.0)THEN
1529 !!$ SurfaceElOnProc = .false.
1530 !!$ exit
1531 !!$ endif
1532 !!$ enddo
1533  IF(epart(ptrtri3%ElemData(4)).EQ.iprocs) surfaceelonproc = .true.
1534  IF(surfaceelonproc.eqv..true.) THEN
1535  ! WRITE(4002,'(10i10)') ABS(NodeFlag(ptrtri3%ElemData(1:3))),&
1536  ! ElFlag(ptrtri3%ElemData(4)),2,1 ! last two are not really used
1537 !!$ IF(SurfaceElOnProc) THEN
1538 !!$ write(4003,*) ABS(NodeFlag(ptrtri3%ElemData(1:6)))
1539 !!$ ENDIF
1540  testnumsurfel = testnumsurfel + 1
1541  elconntable(1:3,testnumsurfel) = abs(nodeflag(ptrtri3%ElemData(1:3)))
1542  elflag_list(testnumsurfel) = elflag(ptrtri3%ElemData(4))
1543  ENDIF
1544  ptrtri3 => ptrtri3%next
1545  ENDDO
1546 
1547  IF(testnumsurfel.NE.numeltet2d(2,iprocs))THEN
1548  print*,'ERROR, number of triangles from link list in Mesh2d'
1549  print*,' different then that in read_patran'
1550  print*,testnumsurfel,numeltet2d(2,iprocs)
1551  stop
1552  ENDIF
1553 
1554 ! write_attr = COM_get_function_handle( 'OUT.write_attribute')
1555 ! set_option = COM_get_function_handle( 'OUT.set_option')
1556 
1557 ! IF(iProcs.EQ.1)THEN
1558 
1559 ! CALL COM_call_function( set_option, 2, 'mode', 'w')
1560 
1561 ! ELSE
1562 
1563 ! CALL COM_call_function( set_option, 2, 'mode', 'a')
1564 
1565 ! ENDIF
1566 
1567 ! do not append process rank -----------------
1568 
1569 ! CALL COM_call_function( set_option, 2, 'rankwidth', '0')
1570 
1571 ! write surface window ------------------------
1572 
1573 ! sur_all = Com_get_attribute_handle( surWin//'.all')
1574 ! CALL COM_call_function( write_attr, 4, 'Rocin/SurfMesh_S.hdf', sur_all,&
1575 ! "solid_surf","00.000000")
1576 
1577 
1578 
1579  ELSE IF(meshtype2d.EQ.4)THEN
1580 
1581 
1582  numnpnew = 0
1583  nullify(ptrhex8)
1584 ! renumber elements nodes
1585 
1586  ptrhex8 => surfmesh_hex8_s_head
1587 
1588  DO WHILE(ASSOCIATED(ptrhex8).eqv..true.)
1589 ! mark that changed surface mesh by changing it to a negative
1590  DO i = 1, 4
1591  glbndnum = ptrhex8%ElemData(i)
1592  IF( nodeflag(glbndnum).GT.0) THEN ! if zero means not on this processor, neg means already renumber node
1593  numnpnew = numnpnew + 1
1594  nodeflag(glbndnum) = -nodeflag(glbndnum) ! overwrite nodeflag here but i need it for non S/F interface
1595  ENDIF
1596  ENDDO
1597  ptrhex8 => ptrhex8%next
1598  ENDDO
1599 
1600 
1601  IF(numelhex2d(2,iprocs).NE.0.AND.numnpnew.EQ.0)THEN
1602  print*,'ERROR, found in surface mesh for Solid/Fluid interface'
1603  print*,'Found ',numelhex2d(2,iprocs),'quads but no nodes'
1604  print*,'Stopping'
1605  stop
1606  ENDIF
1607 
1608  IF(numelhex2d(2,iprocs).EQ.0.AND.numnpnew.NE.0)THEN
1609  print*,'Warning: found in sufrace mesh for Solid/Fluid interface'
1610  print*,'Found ',numnpnew,' nodes but no quads'
1611  ENDIF
1612 
1613 ! write the surface mesh files
1614 
1615 ! WRITE(4002,'(2i10)') NumNpNew, NumElhex2D(2,iProcs)
1616 
1617  nodeflag(:)= abs( nodeflag(:)) ! reset
1618 
1619  ALLOCATE(nodeflag2d(1:2,1:numnpnew)) ! for storing what I'm overwriting
1620  numnpnew = 0
1621 ! renumber elements nodes
1622 
1623  ptrhex8 => surfmesh_hex8_s_head
1624  do while(associated(ptrhex8).eqv..true.)
1625 ! mark that changed surface mesh by changing it to a negative
1626  DO i = 1, 4
1627  glbndnum = ptrhex8%ElemData(i)
1628  IF( nodeflag(glbndnum).GT.0) THEN ! if zero means not on this processor
1629  numnpnew = numnpnew + 1
1630 ! WRITE(4002,'(3(x,e16.9),1i10)') coor(1:3,glbNdNum),NodeFlag(glbNdNum)
1631  nodeflag2d(1,numnpnew) = glbndnum
1632  nodeflag2d(2,numnpnew) = nodeflag(glbndnum)
1633  nodeflag(glbndnum) = - numnpnew
1634  endif
1635  enddo
1636  ptrhex8 => ptrhex8%next
1637  ENDDO
1638 
1639  ! MS: Writing only non-vanishing windows
1640  IF (numnpnew.GT.0) THEN
1641  CALL com_new_attribute( surwin//'.nc', 'n', com_double, 3, 'm')
1642  CALL com_set_size( surwin//'.nc', ini, numnpnew )
1643  CALL com_resize_array(surwin//'.nc', ini, meshcoor, 3)
1644 
1645  DO i = 1, numnpnew
1646  meshcoor(1:3,i) = coor(1:3,nodeflag2d(1,i))
1647  END DO
1648 
1649 
1650  CALL com_new_attribute( surwin//'.bv', 'n', com_integer, 1, '')
1651  ! CALL COM_set_size( surWin//'.bv', iNI, NumNpNew )
1652  CALL com_set_array( surwin//'.bv', ini, nodeflag2d(2,1), 2)
1653 
1654 
1655  ! connectivity
1656  testnumsurfel = numelhex2d(2,iprocs)
1657  print*,'Number of Non-Interacting nodes and elements',numnpnew,testnumsurfel
1658 
1659  CALL com_set_size( surwin//'.:q4', ini, testnumsurfel)
1660  CALL com_resize_array( surwin//'.:q4', ini, elconntable, 4)
1661 
1662 
1663  CALL com_new_attribute( surwin//'.bf2c', 'e', com_integer, 1, '')
1664  CALL com_resize_array(surwin//'.bf2c', ini, elflag_list, 1)
1665 
1666  CALL com_new_attribute( surwin//'.bcflag', 'p', com_integer, 1, '')
1667  CALL com_set_size( surwin//'.bcflag', ini, 1)
1668  CALL com_resize_array(surwin//'.bcflag', ini, bcsurfflagtwo, 1)
1669  bcsurfflagtwo(1) = 2
1670  END IF
1671  ! MS : End
1672 
1673 
1674 ! connectivity
1675  nullify(ptrhex8)
1676  testnumsurfel = 0
1677 
1678  ptrhex8 => surfmesh_hex8_s_head
1679  do while(associated(ptrhex8).eqv..true.)
1680  surfaceelonproc = .false.
1681  IF(epart(ptrhex8%ElemData(5)).EQ.iprocs) surfaceelonproc = .true.
1682  IF(surfaceelonproc.eqv..true.)THEN
1683 ! WRITE(4002,'(10i10)') ABS(NodeFlag(ptrhex8%ElemData(1:4))),&
1684 ! ElFlag(ptrhex8%ElemData(5)),1,1 ! last two are not really used
1685  testnumsurfel = testnumsurfel + 1
1686 
1687  elconntable(1:4,testnumsurfel) = abs(nodeflag(ptrhex8%ElemData(1:4)))
1688  ENDIF
1689  ptrhex8 => ptrhex8%next
1690  ENDDO
1691 
1692  IF(testnumsurfel.NE.numelhex2d(2,iprocs))THEN
1693  print*,'Error: Number of Solid Surface hex elements in linked list'
1694  print*,' does not match that given by NumElhex2D(2,iProcs)'
1695  print*,' In linked list =',testnumsurfel
1696  print*,' NumElhex2D(2,iProcs) =',numelhex2d(2,iprocs)
1697  print*,'Stopping'
1698  ENDIF
1699 
1700  nullify(ptrhex8)
1701 
1702  IF(associated(nodeflag2d).eqv..true.) deallocate(nodeflag2d)
1703 
1704  ENDIF
1705 
1706  ! write surface window ------------------------
1707  CALL com_call_function( set_option, 2, 'mode', 'a')
1708 
1709  sur_all = com_get_attribute_handle( surwin//'.all')
1710 
1711  ! MS: Writing only non-vanishing windows
1712  IF (numnpnew.GT.0) THEN
1713  IF (appendtoggle.EQV..true.) THEN
1714  CALL com_call_function( set_option, 2, 'mode', 'a')
1715  ELSE
1716  CALL com_call_function( set_option, 2, 'mode', 'w')
1717  appendtoggle = .true.
1718  END IF
1719  CALL com_call_function( write_attr, 4, 'Rocin/SurfMesh.'//ichr4, sur_all,&
1720  "isolid","00.000000")
1721  END IF
1722  WRITE (*,*) 'iProc = ', iprocs, ' appendToggle = ', appendtoggle
1723  ! MS : END
1724 
1725 
1726  CALL com_delete_window( surwin)
1727 
1728 
1729  IF(associated(nodeflag2d).eqv..true.) deallocate(nodeflag2d)
1730 
1731 
1732 
1733 end subroutine mesh2d
1734 
Rocout::set_option
void set_option(const char *option_name, const char *option_val)
Set an option for Rocout, such as controlling the output format.
Definition: Rocout.C:552
k
j indices k indices k
Definition: Indexing.h:6
linked_list2
Definition: SourceIMP/utilities/RocfracPrep/Linked_List2.f90:53
write_attr
static void write_attr(std::ostream &os, const COM::Attribute *attr, int i)
Definition: RFC_Window_base_IO_tecplot.C:55
i
blockLoc i
Definition: read.cpp:79
mesh2d
subroutine mesh2d(nprocs, iProcs, ichr4)
Definition: SourceIMP/utilities/RocfracPrep/2Dmesh.f90:53
Aff_transformation_repS2::print
virtual std::ostream & print(std::ostream &os) const
Definition: Aff_transformationS2.h:277
meshdata::user_data_type_procnodelist
Definition: SourceIMP/utilities/RocfracPrep/meshdata.f90:299