-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathreadGridTRIANGLE.f90
262 lines (202 loc) · 8.25 KB
/
readGridTRIANGLE.f90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
! ******************************************************************************
!
! $Id: readGridTRIANGLE.f90,v 1.1 2005/03/05 17:25:53 haselbac Exp $
!
! Filename: readGridTRIANGLE.F90
!
! Purpose: Read 2d grid file in TRIANGLE format.
!
! Description: None.
!
! Input: None.
!
! Output: None.
!
! Notes:
! 1. To convert a TRIANGLE grid, three grid files need to be generated by
! TRIANGLE: a .node, a .ele. and a .edge file. The .edge file is needed
! to construct the boundary edge data structure; note that it will only
! be generated by TRIANGLE if the -e switch is specified.
!
! Author: Andreas Haselbacher
!
! Copyright: (c) 2004 by the University of Illinois
!
! RCS Revision history:
!
! $Log: readGridTRIANGLE.f90,v $
! Revision 1.1 2005/03/05 17:25:53 haselbac
! Initial revision
!
! Revision 1.1 2004/01/13 03:45:22 haselbac
! Initial revision
!
! ******************************************************************************
SUBROUTINE readGridTRIANGLE
USE modError
USE modGlobals
USE modGrid
IMPLICIT NONE
! ******************************************************************************
! Declarations and definitions
! ******************************************************************************
! ==============================================================================
! Local variables
! ==============================================================================
INTEGER :: bType,dummyInteger,ib,ic,ie,iFile,iv,nEdges,v1,v2
INTEGER, DIMENSION(:,:), ALLOCATABLE :: boundMap,c2v
CHARACTER :: choice
CHARACTER*(MAX_STRING_LEN) :: dummyString,iFileName,iFileNameBase
! ******************************************************************************
! Start
! ******************************************************************************
WRITE(STDOUT,'(/,1X,A)') 'Enter file name (without .<n>.<ext>):'
READ(STDIN,'(A)') iFileNameBase
WRITE(STDOUT,'(/)')
WRITE(STDOUT,'(1X,A)') 'Reading grid file in TRIANGLE format...'
! ******************************************************************************
! Read .node file
! ******************************************************************************
iFile = FILE_UNIT_GRID_INPUT
iFileName = TRIM(iFileNameBase)//'.1.node'
OPEN(iFile,FILE=iFileName,FORM="FORMATTED",STATUS="OLD",IOSTAT=errorFlag)
IF ( errorFlag /= NO_ERROR ) THEN
CALL errorHandling(FILE_OPEN_ERROR,iFileName,errorFlag)
END IF ! errorFlag
READ(iFile,*) grid%nVert
! ==============================================================================
! Read coordinates
! ==============================================================================
ALLOCATE(grid%xy(2,grid%nVert),STAT=errorFlag)
IF ( errorFlag /= NO_ERROR ) THEN
CALL errorHandling(ALLOCATE_ERROR,'grid%xy',errorFlag)
END IF ! errorFlag
WRITE(STDOUT,'(3X,A)') 'Coordinates...'
DO iv = 1,grid%nVert
READ(iFile,*) dummyInteger,grid%xy(1,iv),grid%xy(2,iv)
END DO ! iv
CLOSE(iFile,IOSTAT=errorFlag)
IF ( errorFlag /= NO_ERROR ) THEN
CALL errorHandling(FILE_CLOSE_ERROR,iFileName,errorFlag)
END IF ! errorFlag
! ******************************************************************************
! Read .ele file
! ******************************************************************************
iFile = FILE_UNIT_GRID_INPUT
iFileName = TRIM(iFileNameBase)//'.1.ele'
OPEN(iFile,FILE=iFileName,FORM="FORMATTED",STATUS="OLD",IOSTAT=errorFlag)
IF ( errorFlag /= NO_ERROR ) THEN
CALL errorHandling(FILE_OPEN_ERROR,iFileName,errorFlag)
END IF ! errorFlag
READ(iFile,*) grid%nTris
ALLOCATE(grid%tri2v(3,grid%nTris),STAT=errorFlag)
IF ( errorFlag /= NO_ERROR ) THEN
CALL errorHandling(ALLOCATE_ERROR,'grid%tri2v',errorFlag)
END IF ! errorFlag
WRITE(STDOUT,'(3X,A)') 'Connectivity...'
DO ic = 1,grid%nTris
READ(iFile,*) dummyInteger,grid%tri2v(1,ic), &
grid%tri2v(2,ic), &
grid%tri2v(3,ic)
END DO ! ic
grid%nQuads = 0
grid%nCells = grid%nTris
CLOSE(iFile,IOSTAT=errorFlag)
IF ( errorFlag /= NO_ERROR ) THEN
CALL errorHandling(FILE_CLOSE_ERROR,iFileName,errorFlag)
END IF ! errorFlag
! ******************************************************************************
! Read .edge file
! ******************************************************************************
iFile = FILE_UNIT_GRID_INPUT
iFileName = TRIM(iFileNameBase)//'.1.edge'
OPEN(iFile,FILE=iFileName,FORM="FORMATTED",STATUS="OLD",IOSTAT=errorFlag)
IF ( errorFlag /= NO_ERROR ) THEN
CALL errorHandling(FILE_OPEN_ERROR,iFileName,errorFlag)
END IF ! errorFlag
WRITE(STDOUT,'(3X,A)') 'Boundaries...'
! ==============================================================================
! Determine number of boundaries
! ==============================================================================
READ(iFile,*) nEdges
ALLOCATE(boundMap(2,NBOUNDS_MAX),STAT=errorFlag)
IF ( errorFlag /= NO_ERROR ) THEN
CALL errorHandling(ALLOCATE_ERROR,'grid%tri2v',errorFlag)
END IF ! errorFlag
DO ib = 1,NBOUNDS_MAX
boundMap(1,ib) = 0
boundMap(2,ib) = 0
END DO ! ib
DO ie = 1,nEdges
READ(iFile,*) dummyInteger,dummyInteger,dummyInteger,ib
IF ( ib > 0 .AND. ib <= NBOUNDS_MAX ) THEN
boundMap(1,ib) = boundMap(1,ib) + 1
ELSE IF ( ib > NBOUNDS_MAX ) THEN
CALL errorHandling(NBOUNDS_MAX_ERROR)
END IF ! ib
END DO ! ie
grid%nBounds = 0
DO ib = 1,NBOUNDS_MAX
IF ( boundMap(1,ib) /= 0 ) THEN
grid%nBounds = grid%nBounds + 1
boundMap(2,ib) = grid%nBounds
END IF ! boundMap
END DO ! ib
WRITE(STDOUT,'(5X,A,1X,I2,1X,A)') 'Detected',grid%nBounds,'boundaries.'
WRITE(STDOUT,'(5X,A,1X,A)') 'Boundary','nEdges'
DO ib = 1,grid%nBounds
WRITE(STDOUT,'(6X,1X,I2,5X,I4)') ib,boundMap(1,ib)
END DO ! ib
! ==============================================================================
! Allocate memory for boundary type
! ==============================================================================
ALLOCATE(grid%bound(grid%nBounds),STAT=errorFlag)
IF ( errorFlag /= NO_ERROR ) THEN
CALL errorHandling(ALLOCATE_ERROR,'grid%bound',errorFlag)
END IF ! errorFlag
DO ib = 1,grid%nBounds
grid%bound(ib)%nEdges = boundMap(1,ib)
ALLOCATE(grid%bound(ib)%e2v(2,grid%bound(ib)%nEdges),STAT=errorFlag)
IF ( errorFlag /= NO_ERROR ) THEN
CALL errorHandling(ALLOCATE_ERROR,'grid%bound%e2v',errorFlag)
END IF ! errorFlag
grid%bound(ib)%nEdges = 0 ! Reset
END DO ! ib
! ==============================================================================
! Assign boundary conditions to existing boundaries
! ==============================================================================
WRITE(STDOUT,'(5X,A)') 'Assign boundary numbers from following selection:'
WRITE(STDOUT,'(7X,A)') '100 - Inflow'
WRITE(STDOUT,'(7X,A)') '200 - Outflow'
WRITE(STDOUT,'(7X,A)') '300 - No-slip wall'
WRITE(STDOUT,'(7X,A)') '400 - Slip wall'
WRITE(STDOUT,'(7X,A)') '500 - Farfield'
DO ib = 1,grid%nBounds
WRITE(STDOUT,'(5X,A,1X,I2)') 'Enter boundary number for boundary:',ib
READ(STDIN,*) grid%bound(ib)%bType
END DO ! ib
! ==============================================================================
! Rewind file and read again to get boundary edge connectivity
! ==============================================================================
REWIND(iFile)
READ(iFile,*) dummyInteger
DO ie = 1,nEdges
READ(iFile,*) dummyInteger,v1,v2,ib
IF ( ib > 0 ) THEN
grid%bound(ib)%nEdges = grid%bound(ib)%nEdges + 1
grid%bound(ib)%e2v(1,grid%bound(ib)%nEdges) = v1
grid%bound(ib)%e2v(2,grid%bound(ib)%nEdges) = v2
END IF ! ib
END DO ! ie
! ******************************************************************************
! Close file
! ******************************************************************************
CLOSE(iFile,IOSTAT=errorFlag)
IF ( errorFlag /= NO_ERROR ) THEN
CALL errorHandling(FILE_CLOSE_ERROR,iFileName,errorFlag)
END IF ! errorFlag
WRITE(STDOUT,'(1X,A,/)') 'Grid file read successfully.'
! ******************************************************************************
! End
! ******************************************************************************
END SUBROUTINE readGridTRIANGLE