The TELEMAC-MASCARET system  trunk
flusec_gai.f
Go to the documentation of this file.
1 ! *********************
2  SUBROUTINE flusec_gai
3 ! *********************
4 !
5  &(gloseg,dimglo,nseg,npoin,dt,mesh,unsv2d,flodel,flulim,hz,
6  & icla,doplot)
7 !
8 !***********************************************************************
9 ! TELEMAC2D V7P2
10 !***********************************************************************
11 !
14 !
24 !
25 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
39 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
40 !
41  USE bief
44  & fluxlinedata_flusec2,
45  & deja_flusec2,volflux_flusec2,
46  & flux_flusec2,
47  & numberoflines_flusec2,
48  & time_flusec2
49 !
50  USE interface_parallel, ONLY : p_sum
51  IMPLICIT NONE
52 !
53 !!-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
54 !
55  INTEGER, INTENT(IN) :: NSEG,ICLA,NPOIN
56  INTEGER, INTENT(IN) :: DIMGLO
57  INTEGER, INTENT(IN) :: GLOSEG(dimglo,2)
58  DOUBLE PRECISION, INTENT(IN) :: DT
59  TYPE(bief_mesh) :: MESH
60  TYPE(bief_obj), INTENT(IN) :: FLODEL,UNSV2D, HZ
61  DOUBLE PRECISION, INTENT(IN) :: FLULIM(nseg)
62 !
63 !!-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
64 !
65 ! VOLFLUX_FLUSEC2: CUMULATED VOLUME THROUGH SECTIONS
66 ! FLUX_FLUSEC2: FLUX_FLUSEC2 THROUGH CONTROL SECTIONS
67 !
68  INTEGER, PARAMETER :: MAXEDGES=500
69 !
70  INTEGER IERR
71  LOGICAL DOPLOT
72 !
73  INTEGER I,INP
74  INTEGER ISEC
75  INTEGER MYPOS
76  INTEGER MAXNUMBEROFCLASSES
77 !
78  DOUBLE PRECISION, DIMENSION (2) :: SEG1,SEG2
79  DOUBLE PRECISION :: SEGMENTFLUX
80  DOUBLE PRECISION :: SIGN1,SIGN2
81 !
82  DOUBLE PRECISION :: SUMFLUX, SUMVOLFLUX
83 !
84  DOUBLE PRECISION :: SEGXMIN,SEGXMAX
85  DOUBLE PRECISION :: SEGYMIN,SEGYMAX
86  DOUBLE PRECISION,ALLOCATABLE :: FLUXLINES (:,:)
87 !
88 !----------------------------------------------------------------------
89 !
90 ! PART I
91 !
92 ! SEARCH AND SAVE SEGMENTS (FIRST RUN ONLY)
93 !
94 !----------------------------------------------------------------------
95 !
96  IF(.NOT.deja_flusec2) THEN
97 !
98  inp=gai_files(gaiflx)%LU
99  maxnumberofclasses = 20
100  time_flusec2 = 0.0d0
101 !
102 !------- OPEN FLUXLINE FILE
103 !
104  READ(inp,*) numberoflines_flusec2
105 ! ALLOCATE THE FLUXLINES
106  IF (.NOT.ALLOCATED(fluxlines)) THEN
107  ALLOCATE (fluxlines(numberoflines_flusec2,9), stat=ierr)
108  IF(ierr.NE.0) THEN
109  WRITE(lu,*)'FLUSEC_GAI: ERROR OF ALLOCATION:',ierr
110  CALL plante(1)
111  stop
112  ENDIF
113  ENDIF
114 ! READ NODES INTO FLUXLINE
115  DO i = 1,numberoflines_flusec2
116  READ(inp,*) ( fluxlines(i,isec), isec=1,9 )
117  ENDDO
118 !
119  WRITE(lu,*) "FLUXLINES FOUND ",numberoflines_flusec2,
120  & "CLASSES",icla
121 !
122 !------- DYNAMIC ALLOCATION OF FLUX_FLUSEC2, VOLFLUX_FLUSEC2,...
123 !
124  ALLOCATE(flux_flusec2(numberoflines_flusec2,
125  & maxnumberofclasses),stat=ierr)
126  ALLOCATE(volflux_flusec2(numberoflines_flusec2,
127  & maxnumberofclasses),stat=ierr)
128  ALLOCATE(fluxlinedata_flusec2(numberoflines_flusec2),stat=ierr)
129  DO i = 1,numberoflines_flusec2
130  ALLOCATE(fluxlinedata_flusec2(i)%SECTIONIDS(maxedges),
131  & stat=ierr)
132  ALLOCATE(fluxlinedata_flusec2(i)%DIRECTION(maxedges),
133  & stat=ierr)
134  ENDDO
135 !
136  IF(ierr.NE.0) THEN
137  WRITE(lu,200) ierr
138 200 FORMAT(1x,'FLUSEC_GAI: ERROR DURING ALLOCATION
139  & OF MEMORY: ',/,1x,'ERROR CODE: ',1i6)
140  ENDIF
141 !
142 !------ CLEANUP
143 !
144  DO isec =1,numberoflines_flusec2
145  DO i = 1,maxnumberofclasses
146  volflux_flusec2(isec,i) = 0.0d0
147  ENDDO
148  fluxlinedata_flusec2(isec)%NOFSECTIONS = 0
149  ENDDO
150 !
151 !-------LOOP OVER ALL MESH SEGMENTS TO STORE THEM FOR EACH FLUXLINE
152 !
153  DO i = 1,mesh%NSEG
154 !
155  seg1(1) = mesh%X%R(gloseg(i,1))
156  seg1(2) = mesh%Y%R(gloseg(i,1))
157  seg2(1) = mesh%X%R(gloseg(i,2))
158  seg2(2) = mesh%Y%R(gloseg(i,2))
159 ! LOOP OVER ALL FLUXLINES
160  DO isec =1,numberoflines_flusec2
161 !
162 !----------------------------------------------------------
163 !
164 ! SIGN IS USED TO LOOK ON WHICH SIDE OF THE LINE A NODE IS
165 !
166 ! - SIGN IS NEGATIVE IF WE ARE ON THE RIGHT SIDE
167 ! - SIGN IS POSITIVE IF WE ARE ON THE LEFT SIDE
168 ! - SIGN IS ZERO IF WE ARE ON A POINT
169 !
170 !---------------------------------------------------------
171 !
172  sign1 = (seg1(1) - fluxlines(isec,3))*
173  & (fluxlines(isec,2) - fluxlines(isec,4)) -
174  & (seg1(2) - fluxlines(isec,4)) *
175  & (fluxlines(isec,1) - fluxlines(isec,3))
176 
177  sign2 = (seg2(1) - fluxlines(isec,3))*
178  & (fluxlines(isec,2) - fluxlines(isec,4)) -
179  & (seg2(2) - fluxlines(isec,4)) *
180  & (fluxlines(isec,1) - fluxlines(isec,3))
181 !
182 !---------------------------------------------------------
183 !
184 ! THE FLUXLINE SHOULD NEVER CROSS A NODE (BE ZERO)
185 ! IF THIS HAPPENS WE SHIFT THE NODE (RIGHT AND UPWARDS)
186 !
187 !---------------------------------------------------------
188 !
189  IF(sign1.EQ.0.d0) THEN
190  sign1 = (seg1(1)+0.001d0 - fluxlines(isec,3)) *
191  & (fluxlines(isec,2) - fluxlines(isec,4))-
192  & (seg1(2)+0.001d0 - fluxlines(isec,4)) *
193  & (fluxlines(isec,1) - fluxlines(isec,3))
194  ENDIF
195 !
196  IF(sign2.EQ.0.d0) THEN
197  sign2 = (seg2(1)+0.001d0 - fluxlines(isec,3)) *
198  & (fluxlines(isec,2) - fluxlines(isec,4))-
199  & (seg2(2)+0.001d0 - fluxlines(isec,4)) *
200  & (fluxlines(isec,1) - fluxlines(isec,3))
201  ENDIF
202 ! ADD THE SEGMENT ID TO THE NODES
203  IF(sign1*sign2.LT.0.d0) THEN
204 
205  segxmin = min(seg1(1),seg2(1))
206  segxmax = max(seg1(1),seg2(1))
207  segymin = min(seg1(2),seg2(2))
208  segymax = max(seg1(2),seg2(2))
209 !
210  IF((segxmin > fluxlines(isec,5).AND.(segxmax <
211  & fluxlines(isec,7))).AND.
212  & (segymin > fluxlines(isec,6)).AND.(segymax <
213  & fluxlines(isec,8))) THEN
214 !
215  mypos = fluxlinedata_flusec2(isec)%NOFSECTIONS + 1
216  IF(mypos.EQ.maxedges) THEN
217  WRITE(lu,53)
218 53 FORMAT(/,1x,'GAIA IS STOPPED : ',/
219  & ,1x,' REACHED MAXIMUM LIMIT OF EDGES')
220  CALL plante(1)
221  stop
222  ENDIF
223 !
224  fluxlinedata_flusec2(isec)%SECTIONIDS(mypos) = i
225  IF(sign1.GT.0.d0) THEN
226  fluxlinedata_flusec2(isec)%DIRECTION(mypos) = 1
227  ELSE
228  fluxlinedata_flusec2(isec)%DIRECTION(mypos) = -1
229  ENDIF
230  fluxlinedata_flusec2(isec)%NOFSECTIONS = mypos
231 !
232 ! FOR DEBUGGING
233 !
234 ! WRITE(LU,*) 'ADDED SEGMENTS ',
235 ! & I,GLOSEG(I,1),GLOSEG(I,2)
236 ! WRITE(LU,*) 'AT COORDINATES ',
237 ! & SEG1(1),SEG1(2),SEG2(1),SEG2(2)
238 ! WRITE(LU,*) 'SECTIONS FOUND ',
239 ! & FLUXLINEDATA_FLUSEC2(ISEC)%NOFSECTIONS
240  ENDIF
241  ENDIF
242  ENDDO
243  ENDDO
244  ENDIF
245 ! END SEARCH SEGEMENT (DEJA_FLUSEC2)
246  deja_flusec2 = .true.
247 !
248 !----------------------------------------------------------------------
249 !
250 ! PART II
251 !
252 ! ADD THE FLUXES (FLODEL FROM POSITIVE DEPTHS) FOR SEGMENTS
253 !
255 ! CASE! IF A SEGMENT IS SHARED WE NEED THE HALF FLUX_FLUSEC2?
256 !
257 !----------------------------------------------------------------------
258 !
259 ! ONLY INCREASE THE TIME_FLUSEC2 SORRY, THIS INCLUDES THE MORPHOLOGIC FACTOR!
260  IF(icla.EQ.1)THEN
261  time_flusec2 = time_flusec2 + dt
262  ENDIF
263 ! LOOP OVER ALL FLUXLINES
264  DO isec =1,numberoflines_flusec2
265 ! ICLA ARE THE CLASSES
266  flux_flusec2(isec,icla) = 0.0d0
267 ! LOOP OVER SEGMENT
268  DO i = 1,fluxlinedata_flusec2(isec)%NOFSECTIONS
269  segmentflux = fluxlinedata_flusec2(isec)%DIRECTION(i) *
270  & flodel%R(fluxlinedata_flusec2(isec)%SECTIONIDS(i))
271  flux_flusec2(isec,icla) = flux_flusec2(isec,icla)
272  & + segmentflux
273  volflux_flusec2(isec,icla) = volflux_flusec2(isec,icla)
274  & + (segmentflux*dt)
275  ENDDO
276  ENDDO
277 !
278 !----------------------------------------------------------------------
279 !
280 ! PART IIb
281 !
282 ! SCRIPTING AREA ;-)
283 ! ADD A MASSBOX IF YOU LIKE (ADVANCED USERS ONLY) !
284 ! WARNING THIS PART IS NOT PARALLEL
285 !----------------------------------------------------------------------
286 ! IF (DOPLOT) THEN
287 ! MASS = 0.0D0
288 ! MASSTOTAL = 0.0D0
289 ! DO I=1,NPOIN
290 ! POINT(1) = MESH%X%R(I)
291 ! POINT(2) = MESH%Y%R(I)
292 ! IF (POINT(2).GE.31.2) THEN
293 ! MASS = MASS + HZ%R(I) * (1.0D0/ UNSV2D%R(I))
294 ! ENDIF
295 ! MASSTOTAL = MASSTOTAL + HZ%R(I) * (1.0D0/ UNSV2D%R(I))
296 ! ENDDO
297 !
298 ! IF(NCSIZE.GT.1) THEN
299 ! SUMBOX = P_SUM(MASS)
300 ! SUMGLOBAL = P_SUM(MASSTOTAL)
301 ! WRITE(6,FMT=1001) "FLUXLINE_MASSTOTAL",SUMGLOBAL,TIME_FLUSEC2,DT
302 ! WRITE(6,FMT=1001) "FLUXLINE_MASSBOX ",SUMBOX,TIME_FLUSEC2,DT
303 ! ELSE
304 ! WRITE(6,FMT=1001) "FLUXLINE_MASSTOTAL",MASSTOTAL,TIME_FLUSEC2,DT
305 ! WRITE(6,FMT=1001) "FLUXLINE_MASSBOX ",MASS,TIME_FLUSEC2,DT
306 ! ENDIF
307 !
308 !1001 FORMAT (A18,ES22.14,ES22.14,ES22.14)
309 ! ENDIF
310 !----------------------------------------------------------------------
311 !
312 ! PART III
313 !
314 ! SEND THE RESULTS TO FLUXPR_TELEMAC2D
315 !
316 !----------------------------------------------------------------------
317 !
318  IF(doplot) THEN
319  IF(ncsize.GT.1) THEN
320 ! PARALLEL CASE
321 ! PREPARE SINGLE DATA FOR SENDING
322  DO i=1,numberoflines_flusec2
323  sumflux = p_sum(flux_flusec2(i,icla))
324  sumvolflux = p_sum(volflux_flusec2(i,icla))
325  WRITE(lu,fmt=1000) 'FLUXLINE',i,'MYCLASS',icla,
326  & sumflux,sumvolflux,time_flusec2,dt
327  ENDDO
328  ELSE
329 ! SERIAL CASE
330  DO i=1,numberoflines_flusec2
331  WRITE(lu,fmt=1000) 'FLUXLINE',i,'MYCLASS',icla,
332  & flux_flusec2(i,icla),volflux_flusec2(i,icla),time_flusec2,dt
333  ENDDO
334  ENDIF
335  ENDIF
336 !
337 1000 FORMAT (a8,' ',i2,' ',a7,' ',i2,' ',es22.14,es22.14,es22.14,
338  & es22.14)
339 !
340 !----------------------------------------------------------------------
341 !
342  RETURN
343  END
type(bief_file), dimension(maxlu_gai), target gai_files
For storing information on files.
subroutine flusec_gai(GLOSEG, DIMGLO, NSEG, NPOIN, DT, MESH, UNSV2D, FLODEL, FLULIM, HZ, ICLA, DOPLOT)
Definition: flusec_gai.f:8
Definition: bief.f:3