GRIDMAN
grid managment library
vtk.f
Go to the documentation of this file.
1 C> @file convert/vtk.f
2 C> Converting grid into format which can be used for VTK
3 C GRIDMAN, grid managment library. Author: Vladislav Kotov, v.kotov@fz-juelich.de
4 
5 ! Copyright (c) 2017 Forschungszentrum Juelich GmbH
6 ! Vladislav Kotov
7 !
8 ! This file is part of GRIDMAN.
9 !
10 ! GRIDMAN is free software: you can redistribute it and/or modify
11 ! it under the terms of the GNU General Public License as published by
12 ! the Free Software Foundation, either version 3 of the License, or
13 ! (at your option) any later version.
14 !
15 ! GRIDMAN is distributed in the hope that it will be useful,
16 ! but WITHOUT ANY WARRANTY; without even the implied warranty of
17 ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18 ! GNU General Public License for more details.
19 !
20 ! You should have received a copy of the GNU General Public License
21 ! along with GRIDMAN. If not, see <http://www.gnu.org/licenses/>.
22 
23 C>
24 C> \todo VTK output of edges
25 C>
26 
27 C*********************************************************************************
28 C> Convert 2D grid into VTK format
29 C>
30 C> VTK representation includes cells
31 C> and edges which do not belong to any cell
32 C*********************************************************************************
33  SUBROUTINE gridman_2dgrid2vtk(GRID,NPOINTS,NCELLS,NVERTEX,X,
34  s ifirst,lvertex,ivertex,ctype,ierr,
35  s ntheta,theta)
36  USE gridman,ONLY:gridman_grid,gridman_indlist,
37  u gridman_sp,gridman_dp,
38  u gridman_dbg,gridman_unit,gridman_check
39  USE gridman_lib,ONLY:gridman_grid2d_check,gridman_grid2d_chains,
40  u gridman_indlist_deallocate
41  IMPLICIT NONE
42 C> Original grid
43  TYPE(gridman_grid) :: GRID
44 C> Number of points
45  INTEGER(GRIDMAN_SP),INTENT(OUT) :: NPOINTS
46 C> Number of cells
47  INTEGER(GRIDMAN_SP),INTENT(OUT) :: NCELLS
48 C> Total number of vertices
49  INTEGER(GRIDMAN_SP),INTENT(OUT) :: NVERTEX
50 C> Array of coordinates, X(ND,NPOINTS)
51  REAL(GRIDMAN_DP),ALLOCATABLE,INTENT(OUT) :: X(:,:)
52 C> First index which corresponds to each cell in array IVERTEX, IFIRST(NCELLS)
53  INTEGER(GRIDMAN_SP),ALLOCATABLE,INTENT(OUT) :: IFIRST(:)
54 C> Number of vertices of each cell, LVERTEX(NCELLS)
55  INTEGER(GRIDMAN_SP),ALLOCATABLE,INTENT(OUT) :: LVERTEX(:)
56 C> List of points (vertices) belonging to each cell, IVERTEX(NVERTEX)
57  INTEGER(GRIDMAN_SP),ALLOCATABLE,INTENT(OUT) :: IVERTEX(:)
58 C> VTK cell type index, CTYPE(NCELLS)
59  INTEGER,ALLOCATABLE,INTENT(OUT) :: CTYPE(:)
60 C> Error code
61  INTEGER,INTENT(OUT) :: IERR
62 C> Number of nodes in THETA grid
63  INTEGER,INTENT(IN) :: NTHETA
64 C> Cylindrical angle, THETA(NTHETA), [radian]
65 C>
66 C> If NTHETA and THETA(NTHETA) are defined, then axi-symmetric grid is generated
67  REAL(GRIDMAN_DP),INTENT(IN),OPTIONAL :: THETA(ntheta)
68 
69  INTEGER(GRIDMAN_SP) :: IEDGE,IV
70  INTEGER :: ND,IS,RES,IERR0
71  TYPE(gridman_indlist) :: CHAINS
72 
73  IF(gridman_dbg)
74  w WRITE(gridman_unit,*) "Starting GRIDMAN_2DGRID2VTK"
75 
76  ierr=0
77 
78  IF(grid%TYPE.NE.2.OR.grid%PDIM.NE.2.OR.grid%EDIM.NE.2) THEN
79  ierr=100
80  WRITE(gridman_unit,*) "ERROR in GRIDMAN_2DGRID2VTK: ",
81  w "the grid is not of 2D type"
82  WRITE(gridman_unit,*) " TYPE, PDIM, EDIM ",
83  w grid%TYPE,grid%PDIM,grid%EDIM
84  RETURN
85  END IF
86 
87  IF(gridman_check) THEN
88  CALL gridman_grid2d_check(grid,res,ierr)
89  IF(res.NE.0.OR.ierr.GT.0) THEN
90  WRITE(gridman_unit,*) "ERROR in GRIDMAN_2DGRID2VTK: ",
91  w "incorrect GRID object"
92  ierr=100
93  RETURN
94  END IF
95  END IF
96 
97  CALL gridman_grid2d_chains(grid,chains,ierr)
98  IF(ierr.NE.0) THEN
99  WRITE(gridman_unit,*) "GRIDMAN_2DGRID2VTK terminated"
100  RETURN
101  END IF
102 
103  ncells=grid%NCELLS
104  npoints=grid%NPOINTS
105  nvertex=chains%L-chains%N
106 
107  IF(chains%N.NE.ncells.OR.nvertex.LT.3*ncells) THEN
108  WRITE(gridman_unit,*) "ERROR in GRIDMAN_2DGRID2VTK: ",
109  w "internal error after GRIDMAN_GRID2D_CHAINS"
110  WRITE(gridman_unit,*) " L, N, NCELLS ",chains%L,chains%N,ncells
111  ierr=400
112  RETURN
113  END IF
114 
115 C INVERT SEQUENCE OF POINTS IN THE CHAIN
116 C IF SIGNED AREA IS NEGATIVE
117  CALL rearrange_chains(ierr)
118  IF(ierr.NE.0) GOTO 100
119 
120 C ADD EDGES WHICH DO NOT BELONG TO ANY CELL
121  DO iedge=1,grid%NEDGES
122  IF(grid%CELLS(1,iedge).EQ.0.AND.
123  f grid%CELLS(2,iedge).EQ.0) THEN
124  ncells=ncells+1
125  nvertex=nvertex+2
126  END IF
127  END DO
128 
129  IF(ntheta.GT.1.AND.PRESENT(theta)) THEN
130 C INCREASE NUMBER OR CELLS FOR 3D
131 C> If resulting grid is 3D, then edges not belonging to any cell are skipped.
132 C> The axi-symmetric grid is represented in X-Y-Z coordinates
133  ncells=ncells*(ntheta-1)
134  nvertex=2*nvertex*(ntheta-1)
135  npoints=npoints*ntheta
136  nd=3 !3D grid
137  ELSE
138  nd=2 !planar 2D grid
139  END IF
140 
141  ALLOCATE(x(nd,npoints),ifirst(ncells),
142  a lvertex(ncells),ivertex(nvertex),
143  a ctype(ncells),stat=is)
144  IF(is.NE.0) THEN
145  WRITE(gridman_unit,*) "ERROR in GRIDMAN_2DGRID2VTK: ",
146  w "can not perform allocation"
147  WRITE(gridman_unit,*) " NPOINTS, NCELLS, NVERTEX ",
148  w npoints, ncells, nvertex
149  ierr=200
150  RETURN
151  END IF
152 
153  IF(nd.EQ.2) THEN
154  CALL create_2d_grid(ierr)
155  IF(ierr.NE.0) GOTO 100
156  ELSE
157  CALL create_3d_grid(ierr)
158  IF(ierr.NE.0) GOTO 100
159  END IF
160 
161 C SHIFT BECAUSE IN VTK INDICES START FROM 0
162  DO iv=1,nvertex
163  ivertex(iv)=ivertex(iv)-1
164  END DO
165 
166  CALL gridman_indlist_deallocate(chains,ierr)
167 
168  IF(gridman_dbg)
169  w WRITE(gridman_unit,*) "GRIDMAN_2DGRID2VTK completed"
170 
171  RETURN
172 
173  100 CALL gridman_indlist_deallocate(chains,ierr0)
174  CALL deallocate_output
175 
176  CONTAINS
177 
178 C
179  SUBROUTINE create_2d_grid(IERR)
180  INTEGER,INTENT(OUT) :: IERR
181  INTEGER(GRIDMAN_SP) :: ICELL,IEDGE,IV,IE
182 
183  x=grid%X
184  iv=1
185  DO icell=1,grid%NCELLS
186  ifirst(icell)=iv
187  lvertex(icell)=chains%ILAST(icell)-chains%IFIRST(icell)
188  DO ie=chains%IFIRST(icell),chains%ILAST(icell)-1
189  ivertex(iv)=chains%IND(ie)
190  iv=iv+1
191  END DO
192  IF(lvertex(icell).GT.4) THEN
193  ctype(icell)=7 !VTK_POLYGON
194  ELSE IF(lvertex(icell).EQ.4) THEN
195  ctype(icell)=9 !VTK_QUAD
196  ELSE
197  ctype(icell)=5 !VTK_TRIANGLE
198  END IF
199  END DO
200 C ADD EDGES WHICH DO NOT BELONG TO ANY CELL
201  icell=grid%NCELLS
202  DO iedge=1,grid%NEDGES
203  IF(grid%CELLS(1,iedge).EQ.0.AND.
204  f grid%CELLS(2,iedge).EQ.0) THEN
205  icell=icell+1
206  lvertex(icell)=2
207  ctype(icell)=3 !VTK_LINE
208  ifirst(icell)=iv
209  ivertex(iv)=grid%POINTS(1,iedge)
210  iv=iv+1
211  ivertex(iv)=grid%POINTS(2,iedge)
212  iv=iv+1
213  END IF
214  END DO
215  IF(iv-1.NE.nvertex.OR.icell.NE.ncells) THEN
216  WRITE(gridman_unit,*) "ERROR in GRIDMAN_2DGRID2VTK: ",
217  w "internal error"
218  WRITE(gridman_unit,*) "Mismatch of dimensions"
219  WRITE(gridman_unit,*) " NVERTEX, IV ",nvertex,iv-1
220  WRITE(gridman_unit,*) " NCELLS, ICELL ",ncells,icell
221  ierr=400
222  RETURN
223  END IF
224  END SUBROUTINE create_2d_grid
225 
226 C
227  SUBROUTINE create_3d_grid(IERR)
228  INTEGER,INTENT(OUT) :: IERR
229  INTEGER(GRIDMAN_SP) :: IPOINT,ICELL,IE,IV,IEDGE,
230  i ipoint1,icell1,i,npoints0
231  INTEGER :: ITOR
232  REAL(GRIDMAN_DP) :: COSTOR,SINTOR
233  INTRINSIC cos,sin
234 
235  npoints0=grid%NPOINTS
236  ipoint1=0
237  DO itor=1,ntheta
238  costor=cos(theta(itor))
239  sintor=sin(theta(itor))
240  DO ipoint=1,npoints0
241  ipoint1=ipoint1+1
242  x(1,ipoint1)=grid%X(1,ipoint)*costor
243  x(2,ipoint1)=grid%X(1,ipoint)*sintor
244  x(3,ipoint1)=grid%X(2,ipoint)
245  END DO
246  END DO
247  icell1=0
248  iv=1
249  DO itor=1,ntheta-1
250  DO icell=1,grid%NCELLS
251  icell1=icell1+1
252  lvertex(icell1)=chains%ILAST(icell)-chains%IFIRST(icell)
253  ifirst(icell1)=iv
254  IF(lvertex(icell1).EQ.3) THEN
255  ctype(icell1)=13 !VTK_WEDGE
256  ELSE IF(lvertex(icell1).EQ.4) THEN
257  ctype(icell1)=12 !VTK_HEXAHEDRON
258  ELSE
259  WRITE(gridman_unit,*) "ERROR in GRIDMAN_2DGRID2VTK: ",
260  w "unsupported cell type"
261  WRITE(gridman_unit,*)
262  w "3D grids with cells which have ",
263  w "number of corners other than 3 or 4 ",
264  w "are not supported for VTK"
265  WRITE(gridman_unit,*) " ICELL, NCORNERS ",
266  w icell, lvertex(icell1)
267  ierr=100
268  RETURN
269  END IF
270  ie=chains%IFIRST(icell)
271  DO i=0,lvertex(icell1)-1
272  ivertex(iv)=chains%IND(ie+i)+npoints0*(itor-1)
273  iv=iv+1
274  END DO
275  DO i=0,lvertex(icell1)-1
276  ivertex(iv)=chains%IND(ie+i)+npoints0*itor
277  iv=iv+1
278  END DO
279  lvertex(icell1)=lvertex(icell1)*2
280  END DO
281 C ADD EDGES WHICH DO NOT BELONG TO ANY CELL
282  DO iedge=1,grid%NEDGES
283  IF(grid%CELLS(1,iedge).EQ.0.AND.
284  f grid%CELLS(2,iedge).EQ.0) THEN
285  icell1=icell1+1
286  lvertex(icell1)=4
287  ctype(icell1)=9 !VTK_QUAD
288  ifirst(icell1)=iv
289  ivertex(iv)=grid%POINTS(1,iedge)+npoints0*(itor-1)
290  iv=iv+1
291  ivertex(iv)=grid%POINTS(2,iedge)+npoints0*(itor-1)
292  iv=iv+1
293  ivertex(iv)=grid%POINTS(2,iedge)+npoints0*itor
294  iv=iv+1
295  ivertex(iv)=grid%POINTS(1,iedge)+npoints0*itor
296  iv=iv+1
297  END IF
298  END DO
299 
300  END DO
301 C CHECKS
302  IF(iv-1.NE.nvertex.OR.icell1.NE.ncells.OR.
303  f ipoint1.NE.npoints) THEN
304  WRITE(gridman_unit,*) "ERROR in GRIDMAN_2DGRID2VTK: ",
305  w "internal error"
306  WRITE(gridman_unit,*) "Mismatch of dimensions"
307  WRITE(gridman_unit,*) " NVERTEX, IV ",nvertex,iv-1
308  WRITE(gridman_unit,*) " NCELLS, ICELL1 ",ncells,icell1
309  WRITE(gridman_unit,*) " NPOINTS, IPOINT1 ",npoints,ipoint1
310  ierr=400
311  RETURN
312  END IF
313  iv=0
314  DO icell=1,ncells
315  iv=iv+lvertex(icell)
316  END DO
317  IF(iv.NE.nvertex) THEN
318  WRITE(gridman_unit,*) "ERROR in GRIDMAN_2DGRID2VTK: ",
319  w "internal error"
320  WRITE(gridman_unit,*) "Mismatch of dimensions"
321  WRITE(gridman_unit,*) " NVERTEX, IV ",nvertex,iv
322  ierr=400
323  RETURN
324  END IF
325  DO iv=1,nvertex
326  IF(ivertex(iv).GT.npoints) THEN
327  WRITE(gridman_unit,*) "ERROR in GRIDMAN_2DGRID2VTK: ",
328  w "internal error"
329  WRITE(gridman_unit,*) "Wrong point index"
330  WRITE(gridman_unit,*) "IV, IVERTEX, NPOINTS ",
331  w iv, ivertex(iv), npoints
332  END IF
333  END DO
334 
335  END SUBROUTINE create_3d_grid
336 
337 C
338  SUBROUTINE rearrange_chains(IERR)
339  INTEGER,INTENT(OUT) :: IERR
340  INTEGER(GRIDMAN_SP) :: IC,IE,IE1,IE2,IP,N,ITMP,IEE
341  REAL(GRIDMAN_DP) :: X1,Y1,X2,Y2,S
342  INTRINSIC REAL,FLOOR
343 
344  ierr=0
345 
346  DO ic=1,chains%N
347 C CALCULATE SIGNED AREA
348  s=0.
349  ie1=chains%IFIRST(ic)
350  ie2=chains%ILAST(ic)
351  n=ie2-ie1
352  IF(n.LT.3) THEN
353  WRITE(gridman_unit,*) "ERROR in GRIDMAN_2DGRID2VTK: ",
354  w "internal error after GRIDMAN_GRID2D_CHAINS"
355  WRITE(gridman_unit,*) " Wrong number of vertices"
356  WRITE(gridman_unit,*) " ICELL, IFIRST, ILAST ",
357  w ic,chains%IFIRST(ic),chains%ILAST(ic)
358  ierr=400
359  RETURN
360  END IF
361  ip=chains%IND(ie1)
362  x1=grid%X(1,ip)
363  y1=grid%X(2,ip)
364  ie1=ie1+1
365  DO ie=ie1,ie2
366  ip=chains%IND(ie)
367  x2=grid%X(1,ip)
368  y2=grid%X(2,ip)
369  s=s+x1*y2-x2*y1
370  x1=x2
371  y1=y2
372  END DO
373  ie2=ie2-1
374  IF(s.LT.0.) THEN
375 C FLIP THE INDEX ARRAY
376  n=floor(REAL((ie2-ie1+1)/2.))-1
377  iee=ie2
378  DO ie=ie1,ie1+n
379  itmp=chains%IND(ie)
380  chains%IND(ie)=chains%IND(iee)
381  chains%IND(iee)=itmp
382  iee=iee-1
383  END DO
384  END IF
385  END DO
386 
387  END SUBROUTINE rearrange_chains
388 
389 C
390  SUBROUTINE deallocate_output
391 
392  IF(ALLOCATED(x)) THEN
393  DEALLOCATE(x,stat=is)
394  IF(is.NE.0) THEN
395  WRITE(gridman_unit,*) "WARNING from GRIDMAN_2DGRID2VTK: ",
396  w "can not deallocate X"
397  END IF
398  END IF
399  IF(ALLOCATED(ifirst)) THEN
400  DEALLOCATE(ifirst,stat=is)
401  IF(is.NE.0) THEN
402  WRITE(gridman_unit,*) "WARNING from GRIDMAN_2DGRID2VTK: ",
403  w "can not deallocate IFIRST"
404  END IF
405  END IF
406  IF(ALLOCATED(lvertex)) THEN
407  DEALLOCATE(lvertex,stat=is)
408  IF(is.NE.0) THEN
409  WRITE(gridman_unit,*) "WARNING from GRIDMAN_2DGRID2VTK: ",
410  w "can not deallocate LVERTEX"
411  END IF
412  END IF
413  IF(ALLOCATED(ivertex)) THEN
414  DEALLOCATE(ivertex,stat=is)
415  IF(is.NE.0) THEN
416  WRITE(gridman_unit,*) "WARNING from GRIDMAN_2DGRID2VTK: ",
417  w "can not deallocate IVERTEX"
418  END IF
419  END IF
420  IF(ALLOCATED(ctype)) THEN
421  DEALLOCATE(ctype,stat=is)
422  IF(is.NE.0) THEN
423  WRITE(gridman_unit,*) "WARNING from GRIDMAN_2DGRID2VTK: ",
424  w "can not deallocate CTYPE"
425  END IF
426  END IF
427  END SUBROUTINE deallocate_output
428 
429  END SUBROUTINE gridman_2dgrid2vtk
430 
431 
432 C*********************************************************************************
433 C> Translate cell and edge indices in VTK input in the sense of
434 C> \ref GRIDMAN_2DGRID2VTK
435 C*********************************************************************************
436  SUBROUTINE gridman_2dgrid2vtk_index(GRID,N,M,RINDEXES,
437  s sindexes,ierr)
438  USE gridman,ONLY:gridman_grid,gridman_sp,gridman_dp,
439  u gridman_dbg,gridman_unit,gridman_check,
440  u gridman_length
441  USE gridman_lib,ONLY:gridman_grid2d_check
442  IMPLICIT NONE
443  INTRINSIC sum,trim,adjustl
444 C> Original grid
445  TYPE(gridman_grid) :: GRID
446 C> Number of elements in cell scalars array
447  INTEGER(GRIDMAN_SP),INTENT(OUT) :: N
448 C> Number of cell scalars arrays
449  INTEGER,INTENT(OUT) :: M
450 C> Cells scalars array, RINDEXES(N,M)
451 C>
452 C> This array contains indices of cells and of edge not belonging to any cell.
453 C> Zero is used for elements for which indices are not defined.
454 C>
455  REAL(GRIDMAN_DP),ALLOCATABLE,INTENT(OUT) :: RINDEXES(:,:)
456 C> Name of each RINDEXES(:,I), SINDEXES(M)
457 C>
458 C> I<number of index>_<number of column>_<INDEX%COLUMNS>
459 C> Designation of edge indices has same format but starts from E instead of I
460  CHARACTER(GRIDMAN_LENGTH),ALLOCATABLE,INTENT(OUT) :: SINDEXES(:)
461 C> Error code
462  INTEGER,INTENT(OUT) :: IERR
463 
464  INTEGER :: RES,I,J,IS,IR,IR1,IR2
465  INTEGER(GRIDMAN_SP) :: IE,ICELL,IEDGE,NEDGES
466  INTEGER(GRIDMAN_SP),ALLOCATABLE :: IEDGES(:)
467  CHARACTER(4) :: STR
468 
469  IF(gridman_dbg)
470  w WRITE(gridman_unit,*) "Starting GRIDMAN_2DGRID2VTK_INDEX"
471 
472  ierr=0
473 
474  IF(grid%TYPE.NE.2.OR.grid%PDIM.NE.2.OR.grid%EDIM.NE.2) THEN
475  ierr=100
476  WRITE(gridman_unit,*) "ERROR in GRIDMAN_2DGRID2VTK_INDEX: ",
477  w "the grid is not of 2D type"
478  WRITE(gridman_unit,*) " TYPE, PDIM, EDIM ",
479  w grid%TYPE,grid%PDIM,grid%EDIM
480  RETURN
481  END IF
482 
483  IF(gridman_check) THEN
484  CALL gridman_grid2d_check(grid,res,ierr)
485  IF(res.NE.0.OR.ierr.GT.0) THEN
486  WRITE(gridman_unit,*) "ERROR in GRIDMAN_2DGRID2VTK_INDEX: ",
487  w "incorrect GRID object"
488  ierr=100
489  RETURN
490  END IF
491  END IF
492 
493  CALL gridman_grid_count(grid,n,ierr)
494  IF(ierr.NE.0) GOTO 100
495 
496  m=0
497  IF(grid%NCELLINDEX.GT.0) m=sum(grid%CELLINDEX(:)%NINDEX)
498  IF(n.GT.grid%NCELLS.AND.grid%NEDGEINDEX.GT.0) THEN
499  m=m+sum(grid%EDGEINDEX(:)%NINDEX)
500  END IF
501 
502  IF(m.LT.1) THEN
503  m=0
504  n=0
505  RETURN
506  END IF
507 
508  ALLOCATE(rindexes(n,m),sindexes(m),stat=is)
509  IF(is.NE.0) GOTO 200
510 
511  ir=0
512  rindexes=0
513  DO i=1,grid%NCELLINDEX
514  ir1=ir+1
515  DO j=1,grid%CELLINDEX(i)%NINDEX
516  ir=ir+1
517  WRITE(str,'(I4)') i
518  sindexes(ir)='I'//trim(adjustl(str))
519  WRITE(str,'(I4)') j
520  sindexes(ir)=trim(sindexes(ir))//'_'//trim(adjustl(str))
521  sindexes(ir)=trim(sindexes(ir))//'_'//
522  / grid%CELLINDEX(i)%COLUMNS(j)
523  END DO
524  ir2=ir
525  DO ie=1,grid%CELLINDEX(i)%NELEMENTS
526  icell=grid%CELLINDEX(i)%INDEXES(0,ie)
527  IF(icell.GT.0.AND.icell.LE.grid%NCELLS) THEN
528  rindexes(icell,ir1:ir2)=grid%CELLINDEX(i)%INDEXES(1:,ie)
529  END IF
530  END DO
531  END DO
532 
533 C ADD INDICES OF EDGES WHICH DO NOT BELONG TO ANY CELL
534  IF(n.GT.grid%NCELLS.AND.grid%NEDGEINDEX.GT.0) THEN
535  ALLOCATE(iedges(grid%NEDGES),stat=is)
536  IF(is.NE.0) GOTO 210
537  nedges=0
538  iedges=0
539  DO iedge=1,grid%NEDGES
540  IF(grid%CELLS(1,iedge).EQ.0.AND.grid%CELLS(2,iedge).EQ.0) THEN
541  nedges=nedges+1
542  iedges(iedge)=nedges
543  END IF
544  END DO
545  IF(nedges.NE.n-grid%NCELLS) GOTO 400
546 
547  DO i=1,grid%NEDGEINDEX
548  ir1=ir+1
549  DO j=1,grid%EDGEINDEX(i)%NINDEX
550  ir=ir+1
551  WRITE(str,'(I4)') i
552  sindexes(ir)='E'//trim(adjustl(str))
553  WRITE(str,'(I4)') j
554  sindexes(ir)=trim(sindexes(ir))//'_'//trim(adjustl(str))
555  sindexes(ir)=trim(sindexes(ir))//'_'//
556  / grid%EDGEINDEX(i)%COLUMNS(j)
557  END DO
558  ir2=ir
559  DO ie=1,grid%EDGEINDEX(i)%NELEMENTS
560  iedge=grid%EDGEINDEX(i)%INDEXES(0,ie)
561  IF(iedge.GT.0.AND.iedge.LE.grid%NEDGES) THEN
562  iedge=iedges(iedge)
563  IF(iedge.GT.0) THEN
564  rindexes(iedge+grid%NCELLS,ir1:ir2)=
565  = grid%EDGEINDEX(i)%INDEXES(1:,ie)
566  END IF
567  END IF
568  END DO
569  END DO
570  DEALLOCATE(iedges,stat=is)
571  END IF !IF(N.GR.GRID%NCELLS.AND.GRID%NEDGEINDEX.GT.0)
572 
573  IF(gridman_dbg)
574  w WRITE(gridman_unit,*) "GRIDMAN_2DGRID2VTK_INDEX completed"
575 
576  RETURN
577 
578  100 WRITE(gridman_unit,*) "GRIDMAN_2DGRID2VTK_INDEX terminated"
579  RETURN
580  200 WRITE(gridman_unit,*) "ERROR in GRIDMAN_2DGRID2VTK_INDEX: ",
581  w "cannot perform allocation"
582  WRITE(gridman_unit,*) " M, N ",m,n
583  ierr=200
584  RETURN
585  210 WRITE(gridman_unit,*) "ERROR in GRIDMAN_2DGRID2VTK_INDEX: ",
586  w "cannot allocate temporary array"
587  WRITE(gridman_unit,*) " NEDGES ",grid%NEDGES
588  ierr=200
589  RETURN
590  400 WRITE(gridman_unit,*) "ERROR in GRIDMAN_2DGRID2VTK_INDEX: ",
591  w "internal error - musmatch between number of edges"
592  WRITE(gridman_unit,*) " NEDGES, NCOUNT-NCELLS ",
593  w nedges, n-grid%NCELLS
594  ierr=400
595  RETURN
596 
597  END SUBROUTINE gridman_2dgrid2vtk_index
598 
599 C*********************************************************************************
600 C> Write 2D grid and data in VTK ASCII legacy format.
601 C>
602 C> Interface to GRIDMAN_VTK_WRITE which calls GRIDMAN_2DGRID2VTK.
603 C*********************************************************************************
604  SUBROUTINE gridman_vtk_grid2d_write(FNAME,HEADER,
605  s grid,nelements,
606  s ncs,ncv,nps,npv,ierr,
607  s cell_scalar,csname,
608  s cell_vector,cvname,
609  s point_scalar,psname,
610  s point_vector,pvname)
611  USE gridman,ONLY:gridman_grid,gridman_sp,gridman_dp,
612  u gridman_dbg,gridman_unit,gridman_length
613  USE gridman_lib,ONLY:gridman_2dgrid2vtk,gridman_2dgrid2vtk_index,
614  u gridman_grid_count,gridman_vtk_write
615  IMPLICIT NONE
616 C> Name of the output file
617  CHARACTER(*),INTENT(IN) :: FNAME
618 C> String with short description of the data set
619  CHARACTER(*),INTENT(IN) :: HEADER
620 C> 2D grid object
621  TYPE(gridman_grid) :: GRID
622 C> First dimension of arrays CELL_SCALAR, CELL_VECTOR
623 C>
624 C> This quantity must be equal to GRID%NCELLS plus the
625 C> number of "free" edges which do not belong to any cell.
626 C> Calculated by GRIDMAN_GRID_COUNT.
627  INTEGER,INTENT(IN) :: NELEMENTS
628 C> Number of scalar cell data sets
629  INTEGER,INTENT(IN) :: NCS
630 C> Number of vector cell data sets
631  INTEGER,INTENT(IN) :: NCV
632 C> Number of scalar point data sets
633  INTEGER,INTENT(IN) :: NPS
634 C> Number of vector point data sets
635  INTEGER,INTENT(IN) :: NPV
636 C> Error code
637  INTEGER,INTENT(OUT) :: IERR
638 C> Scalar cell data
639  REAL(GRIDMAN_DP),INTENT(IN),OPTIONAL ::
640  r cell_scalar(nelements,ncs)
641 C> Names of the scalar cell data sets
642  CHARACTER(*),INTENT(IN),OPTIONAL :: CSNAME(ncs)
643 C> Vector cell data
644  REAL(GRIDMAN_DP),INTENT(IN),OPTIONAL ::
645  r cell_vector(2,nelements,ncv)
646 C> Names of the vector cell data sets
647  CHARACTER(*),INTENT(IN),OPTIONAL :: CVNAME(ncv)
648 C> Scalar point data
649  REAL(GRIDMAN_DP),INTENT(IN),OPTIONAL ::
650  r point_scalar(grid%NPOINTS,nps)
651 C> Names of the scalar point data sets
652  CHARACTER(*),INTENT(IN),OPTIONAL :: PSNAME(nps)
653 C> Vector cell data
654  REAL(GRIDMAN_DP),INTENT(IN),OPTIONAL ::
655  r point_vector(2,grid%NPOINTS,npv)
656 C> Names of the vector point data sets
657  CHARACTER(*),INTENT(IN),OPTIONAL :: PVNAME(npv)
658 
659  INTEGER(GRIDMAN_SP) :: NPOINTS,NCELLS,NVERTEX,NCOUNT
660  REAL(GRIDMAN_DP),ALLOCATABLE :: X(:,:)
661  INTEGER(GRIDMAN_SP),ALLOCATABLE :: IFIRST(:),
662  i lvertex(:),ivertex(:)
663  INTEGER,ALLOCATABLE :: CTYPE(:)
664 C for indexes
665  INTEGER(GRIDMAN_SP) :: N
666  INTEGER :: M,NCS1,IS,ID
667  REAL(GRIDMAN_DP),ALLOCATABLE :: RINDEXES(:,:)
668  CHARACTER(GRIDMAN_LENGTH),ALLOCATABLE :: SINDEXES(:),
669  c local_csname(:)
670  REAL(GRIDMAN_DP),ALLOCATABLE :: LOCAL_CELL_SCALAR(:,:)
671 
672  IF(gridman_dbg)
673  w WRITE(gridman_unit,*) "Starting GRIDMAN_VTK_GRID2D_WRITE"
674 
675  ierr=0
676 
677  CALL gridman_grid_count(grid,ncount,ierr)
678  IF(ierr.NE.0) GOTO 100
679 
680  IF(PRESENT(cell_scalar).OR.PRESENT(cell_vector)) THEN
681  IF(nelements.LT.ncount) THEN
682  WRITE(gridman_unit,*) "ERROR in GRIDMAN_VTK_GRID2D_WRITE: ",
683  w "data arrays are too small"
684  WRITE(gridman_unit,*) " NELEMENTS ",nelements,
685  w " required ",ncount
686  ierr=100
687  RETURN
688  END IF
689  END IF
690 
691  CALL gridman_2dgrid2vtk(grid,npoints,ncells,nvertex,x,
692  s ifirst,lvertex,ivertex,ctype,ierr,0)
693  IF(ierr.NE.0) GOTO 100
694 
695  CALL gridman_2dgrid2vtk_index(grid,n,m,rindexes,sindexes,ierr)
696  IF(ierr.NE.0) GOTO 100
697 
698  IF(PRESENT(cell_scalar)) THEN
699  ncs1=ncs
700  ELSE
701  ncs1=0
702  END IF
703  IF(m.GT.0) THEN
704  IF(n.NE.ncount) GOTO 400
705  ALLOCATE(local_cell_scalar(n,m+ncs1),
706  a local_csname(m+ncs1),stat=is)
707  IF(is.NE.0) GOTO 200
708  IF(m.GT.0) THEN
709  local_cell_scalar(1:n,1:m)=rindexes(1:n,1:m)
710  local_csname(1:m)=sindexes(1:m)
711  END IF
712  IF(ncs1.GT.0) THEN
713  local_cell_scalar(1:n,m+1:m+ncs1)=cell_scalar(1:n,1:ncs1)
714  IF(PRESENT(csname)) THEN
715  local_csname(m+1:m+ncs1)=csname(1:ncs1)
716  ELSE
717  DO id=1,ncs1
718  WRITE(local_csname(m+id),'(A11,I4.4)') "CELL_SCALAR",id
719  END DO
720  END IF
721  END IF
722  CALL gridman_vtk_write(fname,header,npoints,ncells,nvertex,2,
723  s x,ifirst,lvertex,ivertex,ctype,
724  s m+ncs1,ncv,nps,npv,ierr,
725  s local_cell_scalar,local_csname,
726  s cell_vector,cvname,
727  s point_scalar,psname,point_vector,pvname)
728  ELSE !IF(M+NCS1.GT.0)
729  CALL gridman_vtk_write(fname,header,npoints,ncells,nvertex,2,
730  s x,ifirst,lvertex,ivertex,ctype,
731  s ncs,ncv,nps,npv,ierr,
732  s cell_scalar,csname,cell_vector,cvname,
733  s point_scalar,psname,point_vector,pvname)
734  END IF !IF(M+NCS1.GT.0)
735  IF(ierr.NE.0) THEN
736  CALL deallocate_temporary_arrays
737  WRITE(gridman_unit,*) "GRIDMAN_VTK_GRID2D_WRITE terminated"
738  RETURN
739  END IF
740 
741  CALL deallocate_temporary_arrays
742 
743  IF(gridman_dbg)
744  w WRITE(gridman_unit,*) "GRIDMAN_VTK_GRID2D_WRITE finished"
745 
746  RETURN
747 
748  100 WRITE(gridman_unit,*) "GRIDMAN_VTK_GRID2D_WRITE terminated"
749  RETURN
750  200 WRITE(gridman_unit,*) "ERROR in GRIDMAN_VTK_GRID2D_WRITE: ",
751  w "cannot allocate temporary arrays"
752  WRITE(gridman_unit,*) " NELEMENTS,M+NCS1 ",nelements,m+ncs1
753  CALL deallocate_temporary_arrays
754  ierr=200
755  RETURN
756  400 WRITE(gridman_unit,*) "ERROR in GRIDMAN_VTK_GRID2D_WRITE: ",
757  w "internal error"
758  WRITE(gridman_unit,*) "Counts do not match, N, NCOUNT ",n,ncount
759  ierr=400
760  RETURN
761 
762  CONTAINS
763 C
764  SUBROUTINE deallocate_temporary_arrays
765  INTEGER :: IS
766 
767  IF(ALLOCATED(x)) THEN
768  DEALLOCATE(x,stat=is)
769  IF(is.NE.0) THEN
770  WRITE(gridman_unit,*)
771  w "WARNING from GRIDMAN_VTK_GRID2D_WRITE: ",
772  w "can not deallocate X"
773  END IF
774  END IF
775  IF(ALLOCATED(ifirst)) THEN
776  DEALLOCATE(ifirst,stat=is)
777  IF(is.NE.0) THEN
778  WRITE(gridman_unit,*)
779  w "WARNING from GRIDMAN_VTK_GRID2D_WRITE: ",
780  w "can not deallocate IFIRST"
781  END IF
782  END IF
783  IF(ALLOCATED(lvertex)) THEN
784  DEALLOCATE(lvertex,stat=is)
785  IF(is.NE.0) THEN
786  WRITE(gridman_unit,*)
787  w "WARNING from GRIDMAN_VTK_GRID2D_WRITE: ",
788  w "can not deallocate LVERTEX"
789  END IF
790  END IF
791  IF(ALLOCATED(ivertex)) THEN
792  DEALLOCATE(ivertex,stat=is)
793  IF(is.NE.0) THEN
794  WRITE(gridman_unit,*)
795  w "WARNING from GRIDMAN_VTK_GRID2D_WRITE: ",
796  w "can not deallocate IVERTEX"
797  END IF
798  END IF
799  IF(ALLOCATED(ctype)) THEN
800  DEALLOCATE(ctype,stat=is)
801  IF(is.NE.0) THEN
802  WRITE(gridman_unit,*)
803  w "WARNING from GRIDMAN_VTK_GRID2D_WRITE: ",
804  w "can not deallocate CTYPE"
805  END IF
806  END IF
807  IF(ALLOCATED(rindexes)) THEN
808  DEALLOCATE(rindexes,stat=is)
809  IF(is.NE.0) THEN
810  WRITE(gridman_unit,*)
811  w "WARNING from GRIDMAN_VTK_GRID2D_WRITE: ",
812  w "can not deallocate RINDEXES"
813  END IF
814  END IF
815  IF(ALLOCATED(sindexes)) THEN
816  DEALLOCATE(sindexes,stat=is)
817  IF(is.NE.0) THEN
818  WRITE(gridman_unit,*)
819  w "WARNING from GRIDMAN_VTK_GRID2D_WRITE: ",
820  w "can not deallocate SINDEXES"
821  END IF
822  END IF
823  IF(ALLOCATED(local_csname)) THEN
824  DEALLOCATE(local_csname,stat=is)
825  IF(is.NE.0) THEN
826  WRITE(gridman_unit,*)
827  w "WARNING from GRIDMAN_VTK_GRID2D_WRITE: ",
828  w "can not deallocate LOCAL_CSNAME"
829  END IF
830  END IF
831  IF(ALLOCATED(local_cell_scalar)) THEN
832  DEALLOCATE(local_cell_scalar,stat=is)
833  IF(is.NE.0) THEN
834  WRITE(gridman_unit,*)
835  w "WARNING from GRIDMAN_VTK_GRID2D_WRITE: ",
836  w "can not deallocate LOCAL_CELL_SCALAR"
837  END IF
838  END IF
839 
840  END SUBROUTINE deallocate_temporary_arrays
841 
842  END SUBROUTINE gridman_vtk_grid2d_write
843 
844 C*********************************************************************************
845 C> Convert 3D grid into VTK format.
846 C> Only grid w/o cells is implemented at the moment !!!
847 C*********************************************************************************
848  SUBROUTINE gridman_3dgrid2vtk(GRID,NPOINTS,NCELLS,NVERTEX,X,
849  s ifirst,lvertex,ivertex,ctype,ierr)
850  USE gridman,ONLY:gridman_grid,gridman_sp,gridman_dp,
851  u gridman_dbg,gridman_unit
852  IMPLICIT NONE
853 C> Original grid
854  TYPE(gridman_grid) :: GRID
855 C> Number of points
856  INTEGER(GRIDMAN_SP),INTENT(OUT) :: NPOINTS
857 C> Number of cells
858  INTEGER(GRIDMAN_SP),INTENT(OUT) :: NCELLS
859 C> Total number of vertices
860  INTEGER(GRIDMAN_SP),INTENT(OUT) :: NVERTEX
861 C> Array of coordinates, X(3,NPOINTS)
862  REAL(GRIDMAN_DP),ALLOCATABLE,INTENT(OUT) :: X(:,:)
863 C> First index which corresponds to each cell in array IVERTEX, IFIRST(NCELLS)
864  INTEGER(GRIDMAN_SP),ALLOCATABLE,INTENT(OUT) :: IFIRST(:)
865 C> Number of vertices of each cell, LVERTEX(NCELLS)
866  INTEGER(GRIDMAN_SP),ALLOCATABLE,INTENT(OUT) :: LVERTEX(:)
867 C> List of points (vertices) belonging to each cell, IVERTEX(NVERTEX)
868  INTEGER(GRIDMAN_SP),ALLOCATABLE,INTENT(OUT) :: IVERTEX(:)
869 C> VTK cell type index, CTYPE(NCELLS)
870  INTEGER,ALLOCATABLE,INTENT(OUT) :: CTYPE(:)
871 C> Error code
872  INTEGER,INTENT(OUT) :: IERR
873 
874  INTEGER(GRIDMAN_SP) :: IEDGE,IV,I,IP
875  INTEGER :: IS
876 
877  IF(gridman_dbg)
878  w WRITE(gridman_unit,*) "Starting GRIDMAN_3DGRID2VTK"
879 
880  ierr=0
881  IF(grid%TYPE.NE.3.OR.grid%PDIM.NE.3.OR.grid%EDIM.NE.3) THEN
882  ierr=100
883  WRITE(gridman_unit,*) "ERROR in GRIDMAN_3DGRID2VTK: ",
884  w "the grid is not of FULLERGRID type"
885  WRITE(gridman_unit,*) " TYPE ",grid%TYPE
886  RETURN
887  END IF
888 
889  IF(grid%NCELLS.NE.0) THEN
890  ierr=100
891  WRITE(gridman_unit,*) "ERROR in GRIDMAN_3DGRID2VTK: ",
892  w "at the moment only FULLER grids ",
893  w "w/o cells are implemented"
894  RETURN
895  END IF
896 
897  ncells=grid%NEDGES
898  npoints=grid%NPOINTS
899  nvertex=3*ncells
900 
901  ALLOCATE(x(3,npoints),ifirst(ncells),
902  a lvertex(ncells),ivertex(nvertex),
903  a ctype(ncells),stat=is)
904  IF(is.NE.0) THEN
905  WRITE(gridman_unit,*) "ERROR in GRIDMAN_3DGRID2VTK: ",
906  w "can not perform allocation"
907  WRITE(gridman_unit,*) " NPOINTS, NCELLS, NVERTEX ",
908  w npoints, ncells, nvertex
909  ierr=200
910  RETURN
911  END IF
912 
913  ctype=5 !VTK_TRIANGLE
914  x=grid%X
915  lvertex=3
916  iv=0
917  DO iedge=1,ncells
918  ifirst(iedge)=iv+1
919  DO i=1,3
920  iv=iv+1
921  ip=grid%POINTS(i,iedge)
922  ivertex(iv)=ip
923  END DO
924  END DO
925 
926 C SHIFT BECAUSE IN VTK INDICES START FROM 0
927  DO iv=1,nvertex
928  ivertex(iv)=ivertex(iv)-1
929  END DO
930 
931  END SUBROUTINE gridman_3dgrid2vtk
932 
933 C*********************************************************************************
934 C> Translate edge indices in VTK input in the sense of \ref GRIDMAN_3DGRID2VTK
935 C> (in the current implementation of \ref GRIDMAN_3DGRID2VTK
936 C> only grid edges are treated)
937 C*********************************************************************************
938  SUBROUTINE gridman_3dgrid2vtk_index(GRID,N,M,RINDEXES,
939  s sindexes,ierr)
940  USE gridman,ONLY:gridman_grid,gridman_sp,gridman_dp,
941  u gridman_dbg,gridman_unit,gridman_check,
942  u gridman_length
943  USE gridman_lib,ONLY:gridman_grid_check
944  IMPLICIT NONE
945  INTRINSIC sum,trim,adjustl
946 C> Original grid
947  TYPE(gridman_grid) :: GRID
948 C> Number of elements in cell scalars array
949  INTEGER(GRIDMAN_SP),INTENT(OUT) :: N
950 C> Number of cell scalars arrays
951  INTEGER,INTENT(OUT) :: M
952 C> Cells scalars array, RINDEXES(N,M)
953 C>
954 C> This array contains indexes of edges (in current implementation).
955 C> Zero is used for the elements for which no index is defined.
956  REAL(GRIDMAN_DP),ALLOCATABLE,INTENT(OUT) :: RINDEXES(:,:)
957 C> Name of each RINDEXES(:,I), SINDEXES(M)
958 C>
959 C> E<number of index>_<number of column>_<INDEX%COLUMNS>
960  CHARACTER(GRIDMAN_LENGTH),ALLOCATABLE,INTENT(OUT) :: SINDEXES(:)
961 C> Error code
962  INTEGER,INTENT(OUT) :: IERR
963 
964  INTEGER :: RES,I,J,IS,IR,IR1,IR2
965  INTEGER(GRIDMAN_SP) :: IE,IEDGE
966  CHARACTER(4) :: STR
967 
968  IF(gridman_dbg)
969  w WRITE(gridman_unit,*) "Starting GRIDMAN_3DGRID2VTK_INDEX"
970 
971  ierr=0
972 
973  IF(grid%TYPE.NE.3.OR.grid%PDIM.NE.3.OR.grid%EDIM.NE.3) THEN
974  ierr=100
975  WRITE(gridman_unit,*) "ERROR in GRIDMAN_3DGRID2VTK_INDEX: ",
976  w "the grid is not of 3D type"
977  WRITE(gridman_unit,*) " TYPE, PDIM, EDIM ",
978  w grid%TYPE,grid%PDIM,grid%EDIM
979  RETURN
980  END IF
981 
982  IF(gridman_check) THEN
983  CALL gridman_grid_check(grid,res,ierr)
984  IF(res.NE.0.OR.ierr.GT.0) THEN
985  WRITE(gridman_unit,*) "ERROR in GRIDMAN_3DGRID2VTK_INDEX: ",
986  w "incorrect GRID object"
987  ierr=100
988  RETURN
989  END IF
990  END IF
991 
992  m=0
993  IF(grid%NEDGEINDEX.GT.0) m=sum(grid%EDGEINDEX(:)%NINDEX)
994 
995  n=grid%NEDGES
996  ALLOCATE(rindexes(n,m),sindexes(m),stat=is)
997  IF(is.NE.0) GOTO 200
998 
999  rindexes=0
1000  ir=0
1001  DO i=1,grid%NEDGEINDEX
1002  ir1=ir+1
1003  DO j=1,grid%EDGEINDEX(i)%NINDEX
1004  ir=ir+1
1005  WRITE(str,'(I4)') i
1006  sindexes(ir)='E'//trim(adjustl(str))
1007  WRITE(str,'(I4)') j
1008  sindexes(ir)=trim(sindexes(ir))//'_'//trim(adjustl(str))
1009  sindexes(ir)=trim(sindexes(ir))//'_'//
1010  / grid%EDGEINDEX(i)%COLUMNS(j)
1011  END DO
1012  ir2=ir
1013  DO ie=1,grid%EDGEINDEX(i)%NELEMENTS
1014  iedge=grid%EDGEINDEX(i)%INDEXES(0,ie)
1015  IF(iedge.GT.0.AND.iedge.LE.grid%NEDGES) THEN
1016  rindexes(iedge,ir1:ir2)=grid%EDGEINDEX(i)%INDEXES(1:,ie)
1017  END IF
1018  END DO
1019  END DO
1020 
1021  IF(gridman_dbg)
1022  w WRITE(gridman_unit,*) "GRIDMAN_3DGRID2VTK_INDEX completed"
1023 
1024  RETURN
1025 
1026  200 WRITE(gridman_unit,*) "ERROR in GRIDMAN_3DGRID2VTK_INDEX: ",
1027  w "cannot perform allocation"
1028  WRITE(gridman_unit,*) " M, N ",m,n
1029  ierr=200
1030  RETURN
1031 
1032  END SUBROUTINE gridman_3dgrid2vtk_index
1033 
1034 
1035 C*********************************************************************************
1036 C> Write 3D grid and data in VTK ASCII legacy format.
1037 C> Only grid w/o cells is implemented at the moment !!!
1038 C>
1039 C> Interface to GRIDMAN_VTK_WRITE which calls GRIDMAN_3DGRID2VTK
1040 C*********************************************************************************
1041  SUBROUTINE gridman_vtk_grid3d_write(FNAME,HEADER,
1042  s grid,nelements,
1043  s ncs,ncv,nps,npv,ierr,
1044  s cell_scalar,csname,
1045  s cell_vector,cvname,
1046  s point_scalar,psname,
1047  s point_vector,pvname)
1048  USE gridman,ONLY:gridman_grid,gridman_sp,gridman_dp,
1049  u gridman_dbg,gridman_unit,gridman_length
1050  USE gridman_lib,ONLY:gridman_3dgrid2vtk,gridman_vtk_write,
1051  u gridman_3dgrid2vtk_index
1052  IMPLICIT NONE
1053 C> Name of the output file
1054  CHARACTER(*),INTENT(IN) :: FNAME
1055 C> String with short description of the data set
1056  CHARACTER(*),INTENT(IN) :: HEADER
1057 C> Fuller grid object
1058  TYPE(gridman_grid) :: GRID
1059 C> First dimension of arrays CELL_SCALAR, CELL_VECTOR
1060 C> In current implementation must be equal to GRID%NEDGES
1061  INTEGER(GRIDMAN_SP),INTENT(IN) :: NELEMENTS
1062 C> Number of scalar cell data sets
1063  INTEGER,INTENT(IN) :: NCS
1064 C> Number of vector cell data sets
1065  INTEGER,INTENT(IN) :: NCV
1066 C> Number of scalar point data sets
1067  INTEGER,INTENT(IN) :: NPS
1068 C> Number of vector point data sets
1069  INTEGER,INTENT(IN) :: NPV
1070 C> Error code
1071  INTEGER,INTENT(OUT) :: IERR
1072 C> Scalar cell data
1073  REAL(GRIDMAN_DP),INTENT(IN),OPTIONAL ::
1074  r cell_scalar(nelements,ncs)
1075 C> Names of the scalar cell data sets
1076  CHARACTER(*),INTENT(IN),OPTIONAL :: CSNAME(ncs)
1077 C> Vector cell data
1078  REAL(GRIDMAN_DP),INTENT(IN),OPTIONAL ::
1079  r cell_vector(3,nelements,ncv)
1080 C> Names of the vector cell data sets
1081  CHARACTER(*),INTENT(IN),OPTIONAL :: CVNAME(ncv)
1082 C> Scalar point data
1083  REAL(GRIDMAN_DP),INTENT(IN),OPTIONAL ::
1084  r point_scalar(grid%NPOINTS,nps)
1085 C> Names of the scalar point data sets
1086  CHARACTER(*),INTENT(IN),OPTIONAL :: PSNAME(nps)
1087 C> Vector cell data
1088  REAL(GRIDMAN_DP),INTENT(IN),OPTIONAL ::
1089  r point_vector(3,grid%NPOINTS,npv)
1090 C> Names of the vector point data sets
1091  CHARACTER(*),INTENT(IN),OPTIONAL :: PVNAME(npv)
1092 
1093  INTEGER(GRIDMAN_SP) :: NPOINTS,NCELLS,NVERTEX,NCOUNT
1094  REAL(GRIDMAN_DP),ALLOCATABLE :: X(:,:)
1095  INTEGER(GRIDMAN_SP),ALLOCATABLE :: IFIRST(:),
1096  i lvertex(:),ivertex(:)
1097  INTEGER,ALLOCATABLE :: CTYPE(:)
1098 C for indexes
1099  INTEGER(GRIDMAN_SP) :: N
1100  INTEGER :: M,NCS1,IS,ID
1101  REAL(GRIDMAN_DP),ALLOCATABLE :: RINDEXES(:,:)
1102  CHARACTER(GRIDMAN_LENGTH),ALLOCATABLE :: SINDEXES(:),
1103  c local_csname(:)
1104  REAL(GRIDMAN_DP),ALLOCATABLE :: LOCAL_CELL_SCALAR(:,:)
1105 
1106  IF(gridman_dbg)
1107  w WRITE(gridman_unit,*) "Starting GRIDMAN_VTK_GRID3D_WRITE"
1108 
1109  ierr=0
1110 
1111  ncount=grid%NEDGES
1112  IF(PRESENT(cell_scalar).OR.PRESENT(cell_vector)) THEN
1113  IF(nelements.LT.ncount) THEN
1114  WRITE(gridman_unit,*)
1115  w "ERROR in GRIDMAN_VTK_GRID3D_WRITE: ",
1116  w "data arrays are too small"
1117  WRITE(gridman_unit,*) " NELEMENTS ",nelements,
1118  w " required ",ncount
1119  ierr=100
1120  RETURN
1121  END IF
1122  END IF
1123 
1124  CALL gridman_3dgrid2vtk(grid,npoints,ncells,nvertex,x,
1125  s ifirst,lvertex,ivertex,ctype,ierr)
1126  IF(ierr.NE.0) GOTO 100
1127 
1128  CALL gridman_3dgrid2vtk_index(grid,n,m,rindexes,sindexes,ierr)
1129  IF(ierr.NE.0) GOTO 100
1130 
1131  IF(PRESENT(cell_scalar)) THEN
1132  ncs1=ncs
1133  ELSE
1134  ncs1=0
1135  END IF
1136  IF(m.GT.0) THEN
1137  IF(n.NE.ncount) GOTO 400
1138  ALLOCATE(local_cell_scalar(n,m+ncs1),
1139  a local_csname(m+ncs1),stat=is)
1140  IF(is.NE.0) GOTO 200
1141  IF(m.GT.0) THEN
1142  local_cell_scalar(1:n,1:m)=rindexes(1:n,1:m)
1143  local_csname(1:m)=sindexes(1:m)
1144  END IF
1145  IF(ncs1.GT.0) THEN
1146  local_cell_scalar(1:n,m+1:m+ncs1)=cell_scalar(1:n,1:ncs1)
1147  IF(PRESENT(csname)) THEN
1148  local_csname(m+1:m+ncs1)=csname(1:ncs1)
1149  ELSE
1150  DO id=1,ncs1
1151  WRITE(local_csname(m+id),'(A11,I4.4)') "CELL_SCALAR",id
1152  END DO
1153  END IF
1154  END IF
1155  CALL gridman_vtk_write(fname,header,npoints,ncells,nvertex,3,
1156  s x,ifirst,lvertex,ivertex,ctype,
1157  s m+ncs1,ncv,nps,npv,ierr,
1158  s local_cell_scalar,local_csname,
1159  s cell_vector,cvname,
1160  s point_scalar,psname,point_vector,pvname)
1161  ELSE
1162  CALL gridman_vtk_write(fname,header,npoints,ncells,nvertex,3,
1163  s x,ifirst,lvertex,ivertex,ctype,
1164  s ncs,ncv,nps,npv,ierr,
1165  s cell_scalar,csname,cell_vector,cvname,
1166  s point_scalar,psname,point_vector,pvname)
1167  END IF
1168  IF(ierr.NE.0) THEN
1169  CALL deallocate_temporary_arrays
1170  WRITE(gridman_unit,*) " GRIDMAN_VTK_GRID3D_WRITE terminated"
1171  RETURN
1172  END IF
1173 
1174  CALL deallocate_temporary_arrays
1175 
1176  IF(gridman_dbg)
1177  w WRITE(gridman_unit,*) " GRIDMAN_VTK_GRID3D_WRITE finished"
1178 
1179  RETURN
1180 
1181  100 WRITE(gridman_unit,*) " GRIDMAN_VTK_GRID3D_WRITE terminated"
1182  RETURN
1183  200 WRITE(gridman_unit,*) "ERROR in GRIDMAN_VTK_GRID3D_WRITE: ",
1184  w "cannot allocate temporary arrays"
1185  WRITE(gridman_unit,*) " N,M+NCS1 ",n,m+ncs1
1186  CALL deallocate_temporary_arrays
1187  ierr=200
1188  RETURN
1189  400 WRITE(gridman_unit,*) "ERROR in GRIDMAN_VTK_GRID3D_WRITE: ",
1190  w "internal error"
1191  WRITE(gridman_unit,*) "Counts do not match, N, NCOUNT ",n,ncount
1192  ierr=400
1193  RETURN
1194 
1195  CONTAINS
1196 C
1197  SUBROUTINE deallocate_temporary_arrays
1198  INTEGER :: IS
1199 
1200  IF(ALLOCATED(x)) THEN
1201  DEALLOCATE(x,stat=is)
1202  IF(is.NE.0) THEN
1203  WRITE(gridman_unit,*)
1204  w "WARNING from GRIDMAN_VTK_GRID3D_WRITE: ",
1205  w "can not deallocate X"
1206  END IF
1207  END IF
1208  IF(ALLOCATED(ifirst)) THEN
1209  DEALLOCATE(ifirst,stat=is)
1210  IF(is.NE.0) THEN
1211  WRITE(gridman_unit,*)
1212  w "WARNING from GRIDMAN_VTK_GRID3D_WRITE: ",
1213  w "can not deallocate IFIRST"
1214  END IF
1215  END IF
1216  IF(ALLOCATED(lvertex)) THEN
1217  DEALLOCATE(lvertex,stat=is)
1218  IF(is.NE.0) THEN
1219  WRITE(gridman_unit,*)
1220  w "WARNING from GRIDMAN_VTK_GRID3D_WRITE: ",
1221  w "can not deallocate LVERTEX"
1222  END IF
1223  END IF
1224  IF(ALLOCATED(ivertex)) THEN
1225  DEALLOCATE(ivertex,stat=is)
1226  IF(is.NE.0) THEN
1227  WRITE(gridman_unit,*)
1228  w "WARNING from GRIDMAN_VTK_GRID3D_WRITE: ",
1229  w "can not deallocate IVERTEX"
1230  END IF
1231  END IF
1232  IF(ALLOCATED(ctype)) THEN
1233  DEALLOCATE(ctype,stat=is)
1234  IF(is.NE.0) THEN
1235  WRITE(gridman_unit,*)
1236  w "WARNING from GRIDMAN_VTK_GRID3D_WRITE: ",
1237  w "can not deallocate CTYPE"
1238  END IF
1239  END IF
1240  IF(ALLOCATED(rindexes)) THEN
1241  DEALLOCATE(rindexes,stat=is)
1242  IF(is.NE.0) THEN
1243  WRITE(gridman_unit,*)
1244  w "WARNING from GRIDMAN_VTK_GRID3D_WRITE: ",
1245  w "can not deallocate RINDEXES"
1246  END IF
1247  END IF
1248  IF(ALLOCATED(sindexes)) THEN
1249  DEALLOCATE(sindexes,stat=is)
1250  IF(is.NE.0) THEN
1251  WRITE(gridman_unit,*)
1252  w "WARNING from GRIDMAN_VTK_GRID3D_WRITE: ",
1253  w "can not deallocate SINDEXES"
1254  END IF
1255  END IF
1256  IF(ALLOCATED(local_csname)) THEN
1257  DEALLOCATE(local_csname,stat=is)
1258  IF(is.NE.0) THEN
1259  WRITE(gridman_unit,*)
1260  w "WARNING from GRIDMAN_VTK_GRID3D_WRITE: ",
1261  w "can not deallocate LOCAL_CSNAME"
1262  END IF
1263  END IF
1264  IF(ALLOCATED(local_cell_scalar)) THEN
1265  DEALLOCATE(local_cell_scalar,stat=is)
1266  IF(is.NE.0) THEN
1267  WRITE(gridman_unit,*)
1268  w "WARNING from GRIDMAN_VTK_GRID3D_WRITE: ",
1269  w "can not deallocate LOCAL_CELL_SCALAR"
1270  END IF
1271  END IF
1272  END SUBROUTINE deallocate_temporary_arrays
1273 
1274  END SUBROUTINE gridman_vtk_grid3d_write
gridman_grid_check
subroutine gridman_grid_check(GRID, RES, IERR)
Check consistency of the grid data.
Definition: grid1.f:289
gridman::gridman_dp
integer, parameter, public gridman_dp
Kind parameter for real numbers.
Definition: gridman.f:93
gridman_lib
Explicit interfaces to GRIDMAN subroutines and functions.
Definition: gridman.f:251
gridman_grid2d_chains
subroutine gridman_grid2d_chains(GRID, CHAINS, IERR)
Find closed chain of points which form each cell.
Definition: grid2d.f:266
gridman_3dgrid2vtk_index
subroutine gridman_3dgrid2vtk_index(GRID, N, M, RINDEXES, SINDEXES, IERR)
Translate edge indices in VTK input in the sense of GRIDMAN_3DGRID2VTK (in the current implementation...
Definition: vtk.f:940
gridman_grid_count
subroutine gridman_grid_count(GRID, NCOUNT, IERR)
Return the number of cells plus the number of edges not belonging to any cell.
Definition: grid1.f:1360
gridman_grid2d_check
subroutine gridman_grid2d_check(GRID, RES, IERR)
Check correctness of the 2D grid object.
Definition: grid2d.f:37
gridman_indlist_deallocate
subroutine gridman_indlist_deallocate(INDLIST, IERR)
Deallocate list of indices.
Definition: indlist.f:106
gridman::gridman_grid
Data-type which describes a grid as a set of edges, methods in grid.f.
Definition: gridman.f:168
gridman_vtk_grid2d_write
subroutine gridman_vtk_grid2d_write(FNAME, HEADER, GRID, NELEMENTS, NCS, NCV, NPS, NPV, IERR, CELL_SCALAR, CSNAME, CELL_VECTOR, CVNAME, POINT_SCALAR, PSNAME, POINT_VECTOR, PVNAME)
Write 2D grid and data in VTK ASCII legacy format.
Definition: vtk.f:611
gridman_3dgrid2vtk
subroutine gridman_3dgrid2vtk(GRID, NPOINTS, NCELLS, NVERTEX, X, IFIRST, LVERTEX, IVERTEX, CTYPE, IERR)
Convert 3D grid into VTK format. Only grid w/o cells is implemented at the moment !!! ...
Definition: vtk.f:850
gridman::gridman_indlist
Data-type which describes lists of elements with variable number of indices for each element...
Definition: gridman.f:225
gridman_2dgrid2vtk_index
subroutine gridman_2dgrid2vtk_index(GRID, N, M, RINDEXES, SINDEXES, IERR)
Translate cell and edge indices in VTK input in the sense of GRIDMAN_2DGRID2VTK.
Definition: vtk.f:438
gridman
Definition of data types, global constants and variables.
Definition: gridman.f:83
gridman_vtk_grid3d_write
subroutine gridman_vtk_grid3d_write(FNAME, HEADER, GRID, NELEMENTS, NCS, NCV, NPS, NPV, IERR, CELL_SCALAR, CSNAME, CELL_VECTOR, CVNAME, POINT_SCALAR, PSNAME, POINT_VECTOR, PVNAME)
Write 3D grid and data in VTK ASCII legacy format. Only grid w/o cells is implemented at the moment !...
Definition: vtk.f:1048
gridman_vtk_write
subroutine gridman_vtk_write(FNAME, HEADER, NPOINTS, NCELLS, NVERTEX, ND, X, IFIRST, LVERTEX, IVERTEX, CTYPE, NCS, NCV, NPS, NPV, IERR, CELL_SCALAR, CSNAME, CELL_VECTOR, CVNAME, POINT_SCALAR, PSNAME, POINT_VECTOR, PVNAME)
Write grid and data in VTK ASCII legacy format (interface)
Definition: vtk.f:34
gridman_2dgrid2vtk
subroutine gridman_2dgrid2vtk(GRID, NPOINTS, NCELLS, NVERTEX, X, IFIRST, LVERTEX, IVERTEX, CTYPE, IERR, NTHETA, THETA)
Convert 2D grid into VTK format.
Definition: vtk.f:36
gridman::gridman_sp
integer, parameter, public gridman_sp
Kind parameter for integer numbers.
Definition: gridman.f:95
Generated on Thu Oct 5 2017 09:49:56 for GRIDMAN by   doxygen 1.8.9.1