UPP  001
 All Data Structures Files Functions Pages
SMOOTH.f
Go to the documentation of this file.
1 
5 
32 !**********************************************************************
33 !**********************************************************************
34 
35  SUBROUTINE smooth (FIELD,HOLD,IX,IY,SMTH)
36 
37 !**********************************************************************
38 !**********************************************************************
39 
40  implicit none
41 
42  integer :: i1, i2, j, it, i, ix, iy
43  real :: smth1, smth, smth2, smth3, smth4, smth5
44  real :: sum1, sum2
45  REAL field(ix,iy), hold (ix,2)
46  smth1 = 0.25 * smth * smth
47  smth2 = 0.5 * smth * (1.-smth)
48  smth3 = (1.-smth) * (1.-smth)
49  smth4 = (1.-smth)
50  smth5 = 0.5 * smth
51  i1 = 2
52  i2 = 1
53  DO j=2,iy-1
54  it = i1
55  i1 = i2
56  i2 = it
57 !$omp parallel do private(i,sum1,sum2)
58  DO i = 2,ix-1
59  sum1 = field(i-1,j+1) + field(i-1,j-1) &
60  + field(i+1,j+1) + field(i+1,j-1)
61  sum2 = field(i ,j+1) + field(i+1,j ) &
62  + field(i ,j-1) + field(i-1,j )
63  hold(i,i1) = smth1*sum1 + smth2*sum2 + smth3*field(i,j)
64  ENDDO
65  IF (j > 2) then
66 !$omp parallel do private(i)
67  DO i=2,ix-1
68  field(i,j-1) = hold(i,i2)
69  ENDDO
70  endif
71  ENDDO
72 
73 !$omp parallel do private(i)
74  DO i = 2,ix-1
75  field(i,iy-1) = hold(i,i1)
76  ENDDO
77 
78  DO i = 2,ix-1
79  field(i,1) = smth4 * field(i,1) &
80  + smth5 * (field(i-1,1) + field(i+1,1))
81  field(i,iy) = smth4 * field(i,iy) &
82  + smth5 * (field(i-1,iy) + field(i+1,iy))
83  ENDDO
84 
85  DO j = 2,iy-1
86  field(1,j) = smth4 * field(1,j) &
87  + smth5 * (field(1,j-1) + field(1,j+1))
88  field(ix,j) = smth4 * field(ix,j) &
89  + smth5 * (field(ix,j-1) + field(ix,j+1))
90  ENDDO
91 
92  RETURN
93  END
97 
124 !**********************************************************************
125 !**********************************************************************
126 
127  SUBROUTINE smoothc (FIELD,HOLD,IX,IY,SMTH)
128 
129 !**********************************************************************
130 !**********************************************************************
131 
132  implicit none
133 
134  integer :: i1, i2, j, it, i, ix, iy, im1, ip1
135  real :: smth1, smth, smth2, smth3, smth4, smth5
136  real :: sum1, sum2
137  REAL field(ix,iy), hold (ix,2)
138  integer :: iw(ix), ie(ix)
139 !
140  smth1 = 0.25 * smth * smth
141  smth2 = 0.5 * smth * (1.-smth)
142  smth3 = (1.-smth) * (1.-smth)
143  smth4 = (1.-smth)
144  smth5 = 0.5 * smth
145 !
146  do i=2,ix-1
147  ie(i) = i + 1
148  iw(i) = i - 1
149  enddo
150  ie(ix) = 1
151  iw(1) = ix
152 !
153  i1 = 2
154  i2 = 1
155  DO j=2,iy-1
156  it = i1
157  i1 = i2
158  i2 = it
159 !$omp parallel do private(i,sum1,sum2,ip1,im1)
160  DO i = 1,ix
161  ip1 = ie(i)
162  im1 = iw(i)
163  sum1 = field(im1,j+1) + field(im1,j-1) &
164  + field(ip1,j+1) + field(ip1,j-1)
165  sum2 = field(i ,j+1) + field(ip1,j ) &
166  + field(i ,j-1) + field(im1,j )
167  hold(i,i1) = smth1*sum1 + smth2*sum2 + smth3*field(i,j)
168  ENDDO
169  IF (j > 2) then
170 !$omp parallel do private(i)
171  DO i=1,ix
172  field(i,j-1) = hold(i,i2)
173  ENDDO
174  endif
175  ENDDO
176 
177 !$omp parallel do private(i)
178  DO i = 1,ix
179  field(i,iy-1) = hold(i,i1)
180  ENDDO
181 
182  DO i = 1,ix
183  ip1 = ie(i)
184  im1 = iw(i)
185  field(i,1) = smth4 * field(i,1) &
186  + smth5 * (field(im1,1) + field(ip1,1))
187  field(i,iy) = smth4 * field(i,iy) &
188  + smth5 * (field(im1,iy) + field(ip1,iy))
189  ENDDO
190 
191 
192  RETURN
193  END