The TELEMAC-MASCARET system  trunk
cvsp_compress_dp_gaia.f
Go to the documentation of this file.
1 ! ***************************
2  SUBROUTINE cvsp_compress_dp_gaia
3 ! ***************************
4 !
5  &(j, threshold)
6 !
7 !***********************************************************************
8 ! GAIA V8P1 14/03/2013
9 !***********************************************************************
10 !
16 !
20 !
25 !
26 !!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
29 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
30 !
33 !
34  USE cvsp_outputfiles_gaia, ONLY: cp
35 
37  IMPLICIT NONE
38 !
39 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
40 !
41  INTEGER, INTENT(IN) :: J
42  DOUBLE PRECISION, INTENT(IN) :: THRESHOLD
43 !
44 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
45 !
46  INTEGER I,K, MARKERMAX, MARKERCNT, TTT, NNN, JG
47  INTEGER MAXPOS, M, MARKERMAXOLD, MARKERMAXVERYOLD
48  INTEGER MARKER(pro_max_max), MARKERTEMP(pro_max_max)
49 !
50 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
51 !
52 ! USED TO MARK NODES THAT WILL BE KEPT
53 !
54  DOUBLE PRECISION LOSS(pro_max_max)
55 !
56 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
57 ! STORES THE FRACTION ERRORS THAT WILL OCCURE IF THE POINT IS ELEMINATED FROM CURRENT PROFILE
58 !
59  DOUBLE PRECISION MAXDIST
60 !
61 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
62 ! STORES THE MAXIMUM DISTANCE OF ANY NODE IN THE CURRENT LOOP
63 !
64  DOUBLE PRECISION FI, FJ, FK, DI, DJ, DK, THRESH
65 !
66 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
67 !
68 !
69 !-----------------------------------------------------------------------
70 ! PARALLEL: LOCAL TO GLOBAL
71 !-----------------------------------------------------------------------
72 !
73  jg = j
74  IF (ncsize.GT.1) jg = mesh%KNOLG%I(j)
75 !
76 !--------------------------------------------------------------------------
77 ! DEBUG: INIT OUTPUT
78 !-----------------------------------------------------------------------
79 !
80  IF(cvsp_db_gaia(jg,0)) CALL cvsp_p_gaia('./','V_H',jg)
81 
82 !--------------------------------------------------------------------------
83 ! INIT
84 !-----------------------------------------------------------------------
85 !
86  IF(pro_max(j) <= 2) RETURN
87 !
88 !-----------------------------------------------------------------------
89 !FIRST AND LAST POINT WILL ALWAY BE KEPT
90 !-----------------------------------------------------------------------
91 !
92  markermax = 2 !MAXIMUM USED INDEX IN MARKER ARRAY
93  marker(1) = 1 !FIRST WILL ALWAYS BE KEPT
94  marker(2) = pro_max(j) !LAST WILL ALWAYS BE KEPT
95  markercnt = 1 !MAXIMUM USED INDEX IN MARKERTEMP ARRAY
96 !
97  thresh = threshold
98 !
99 !-------------------------------------------------------------------------
100 ! EXTEND THRESHOLD IF NECESSARY
101 !-----------------------------------------------------------------------
102 !
103  DO nnn = 1,1 !4
104  IF (nnn > 1) WRITE(lu,*) 'COMPRESS', j, nnn
105 !
106  thresh = thresh * (10**(1-nnn))
107  markermaxveryold = markermax
108 !
109 !--------------------------------------------------------------------------
110 ! ITERATE UNTIL NOTHING CHANGES ANYMORE
111 !-----------------------------------------------------------------------
112 !
113  DO ttt = 1, pro_max(j) - 2 ! THEROTICAL MAXIMUM NUMBER OF ITERATIONS
114  markermaxold = markermax
115  markercnt = 1
116  markertemp(markercnt) = 1
117 !
118 !--------------------------------------------------------------------------
119 ! LOOP OVER ALL SECTIONS BETWEEN 2 MARKED NODES
120 !-----------------------------------------------------------------------
121 !
122  DO i = 1, markermax-1
123  maxdist = 0 !INITS THE MAXIMUM FRACTION ERROR
124  maxpos = -1 !INITS THE NODE WHICH PRODUCES THE MAXIMUM FRACTION ERROR
125  IF (marker(i+1)-marker(i) >= 2 ) THEN
126 !
127 !--------------------------------------------------------------------------
128 ! LOOP OVER ALL UNMARKED NODES IN BETWEEN 2 MARKED NODES
129 !-----------------------------------------------------------------------
130 !
131  DO m = marker(i) + 1 , marker(i + 1) - 1
132 !
133 !--------------------------------------------------------------------------
134 ! HOW MUCH VOLUME=FRACTION IS LOST IF WE ELIMINATE THIS PROFILEPOINT
135 ! USING "VOLUME" !!! ORIGINAL DOUGLAS-PEUKER: DISTANCE TO INTERCONNECTION !!!
136 ! "ERROR TRIANGLE VOLUME" IS CALCULATED BY GAUSSIAN POLYGON FORMULA!
137 !-----------------------------------------------------------------------
138 !
139  loss(m) = 0.d0
140 !
141  DO k = 1, nsicla
142  IF (nnn.GE.5) THEN
143  fi = pro_f(j,m-1,k)
144  fj = pro_f(j,m+1,k)
145  fk = pro_f(j,m,k)
146  di = pro_d(j,m-1,k)
147  dj = pro_d(j,m+1,k)
148  dk = pro_d(j,m,k)
149  ELSE
150  fi = pro_f(j,marker(i),k)
151  fj = pro_f(j,marker(i+1),k)
152  fk = pro_f(j,m,k)
153  di = pro_d(j,marker(i),k)
154  dj = pro_d(j,marker(i+1),k)
155  dk = pro_d(j,m,k)
156  ENDIF
157 !
158  loss(m) = loss(m) +
159  & abs(0.5d0 * ((fi+fj) * (di-dj) +
160  & (fj+fk) * (dj-dk) +
161  & (fk+fi) * (dk-di)))
162 
163  ENDDO !K
164 !
165  IF(loss(m).GT.maxdist) THEN
166  maxdist = loss(m)
167  maxpos = m
168  ENDIF
169 !
170  ENDDO !M
171 !
172 !-----------------------------------------------------------------------
173 ! IF ANY POINT IS TO FAR OUT OF RANGE: ADD IT TO THE MARKER LIST
174 !-----------------------------------------------------------------------
175 !
176  IF(maxpos > -1 .AND. maxdist > thresh) THEN
177  markercnt = markercnt + 1
178  markertemp(markercnt) = maxpos
179  ENDIF
180 
181  ENDIF
182 !
183 !-----------------------------------------------------------------------
184 ! ADD THE ENDPOINT OF THIS SECTION
185 !-----------------------------------------------------------------------
186 !
187  markercnt = markercnt + 1
188  markertemp(markercnt) = marker(i+1)
189 !
190  ENDDO !I
191 !
192  DO i = 1, markercnt
193  marker(i) = markertemp(i)
194  ENDDO
195  markermax = markercnt
196 !
197  IF (markermax - markermaxold == 0 ) EXIT !STOP ITERATION, AS NOTHING CHANGED!
198  ENDDO ! TTT
199  IF (markermax - markermaxveryold == 0 ) EXIT !STOP ITERATION, AS NOTHING CHANGED!
200  ENDDO ! NNN
201 !
202 !--------------------------------------------------------------------------
203 ! RECREATE THE SORTING PROFILE WITH LESSER NUMBER OF SECTIONS
204 !-----------------------------------------------------------------------
205 !
206  DO k = 1, nsicla
207  DO i = 1, markermax
208  pro_f(j,i,k) = pro_f(j,marker(i),k)
209  pro_d(j,i,k) = pro_d(j,marker(i),k)
210  ENDDO !I
211  ENDDO !K
212 !
213  pro_max(j) = markermax
214 !
215 !--------------------------------------------------------------------------
216 ! BRUTFORCE COMPRESSION IN CASE OF EXCEPTIONAL FRAGMENTATION
217 !-----------------------------------------------------------------------
218 !
219  IF(pro_max(j) > pro_max_max-4*nsicla-4) THEN
220  IF(cp)
221  & WRITE(lu,*) 'CVSP_COMPRESS_DP RESIGNS AND CALLS COMPRESS_BRUT:'
223  ENDIF
224 !
225 !--------------------------------------------------------------------------
226 ! DEBUG: FINAL OUTPUT
227 !-----------------------------------------------------------------------
228 !
229  IF(cvsp_db_gaia(jg,0)) CALL cvsp_p_gaia('./','V_I',jg)
230 !
231 !-----------------------------------------------------------------------
232 !
233  RETURN
234  END SUBROUTINE cvsp_compress_dp_gaia
subroutine cvsp_compress_brut_gaia(J)
subroutine cvsp_p_gaia(PATH_PRE, FILE_PRE, JG)
Definition: cvsp_p_gaia.f:7
logical cp
Logical for debug printouts.
subroutine cvsp_compress_dp_gaia(J, THRESHOLD)
integer, target nsicla
Number of sediment classes of bed material (less than NISCLM)
double precision, dimension(:,:,:), allocatable, target pro_f
Vertical sorting profile: fraction for each layer, class, point.
double precision, dimension(:,:,:), allocatable, target pro_d
Vertical sorting profile: depth for each layer, class, point.
integer pro_max_max
Maximum Number of Profile SECTIONS.
logical function cvsp_db_gaia(J_GLOBAL, TIMESTAMP)
Definition: cvsp_db_gaia.f:7
type(bief_mesh), target mesh
Mesh structure.
integer, dimension(:), allocatable pro_max
Maximum layer number in a vertical sorting profile for each point.