V3FIT
alias.f
1  SUBROUTINE alias_par(gcons, ztemp, gcs, gsc, gcc, gss)
2  USE vmec_main
3  USE realspace, ONLY:ireflect_par, psqrts
4  USE parallel_include_module
5  IMPLICIT NONE
6 C-----------------------------------------------
7 C D u m m y A r g u m e n t s
8 C-----------------------------------------------
9  REAL(dp), PARAMETER :: p5 = 0.5_dp
10  REAL(dp), DIMENSION(nzeta,ntheta3,ns), INTENT(out) :: gcons
11  REAL(dp), DIMENSION(nzeta,ntheta3,ns), INTENT(in) :: ztemp
12  REAL(dp), DIMENSION(0:ntor,0:mpol1,ns) :: gcs, gsc, gcc, gss
13 C-----------------------------------------------
14 C L o c a l V a r i a b l e s
15 C-----------------------------------------------
16  INTEGER :: m, i, ir, jk, jka, n, k, js, l, j
17  INTEGER :: nsmin, nsmax
18  INTEGER :: jcount, kk
19  REAL(dp), DIMENSION(:,:,:), ALLOCATABLE :: work
20  REAL(dp), DIMENSION(:,:,:), ALLOCATABLE :: gcona
21  REAL(dp) :: talon, taloff
22 C-----------------------------------------------
23  CALL second0(talon)
24 
25  nsmin=tlglob; nsmax=t1rglob
26 
27  ALLOCATE (work(4,nzeta,ns), gcona(nzeta,ntheta3,ns))
28 
29  gcons(:,:,nsmin:nsmax) = 0
30  gcona(:,:,nsmin:nsmax) = 0
31 
32  gcs(:ntor,:mpol1,nsmin:nsmax) = 0
33  gsc(:ntor,:mpol1,nsmin:nsmax) = 0
34  gcc(:ntor,:mpol1,nsmin:nsmax) = 0
35  gss(:ntor,:mpol1,nsmin:nsmax) = 0
36 
37  !BEGIN M-LOOP
38  DO js = nsmin, nsmax
39  DO m = 1, mpol1 - 1
40 
41  work(:,:,js) = 0
42  DO i = 1, ntheta2
43  DO k = 1, nzeta
44  work(1,k,js) = work(1,k,js) + ztemp(k,i,js)*cosmui(i,m)
45  work(2,k,js) = work(2,k,js) + ztemp(k,i,js)*sinmui(i,m)
46  END DO
47 
48  IF (.not.lasym) cycle
49  ir = ntheta1 + 2 - i
50  IF (i .eq. 1) ir = 1
51  DO k = 1, nzeta
52  kk=ireflect_par(k)
53  work(3,k,js) = work(3,k,js) +
54  1 ztemp(kk,ir,js)*cosmui(i,m)
55  work(4,k,js) = work(4,k,js) +
56  1 ztemp(kk,ir,js)*sinmui(i,m)
57  END DO
58  END DO
59 
60  IF(js.GT.1) THEN
61  DO n = 0, ntor
62  DO k = 1, nzeta
63  IF (.not.lasym) THEN
64  gcs(n,m,js) = gcs(n,m, js) + tcon(js)*work(1,k,js)*
65  1 sinnv(k,n)
66  gsc(n,m,js) = gsc(n,m,js) + tcon(js)*work(2,k,js)*
67  1 cosnv(k,n)
68  ELSE
69  gcs(n,m,js) = gcs(n,m,js) + p5*tcon(js)*sinnv(k,n)*
70  1 (work(1,k,js)-work(3,k,js))
71  gsc(n,m,js) = gsc(n,m,js) + p5*tcon(js)*cosnv(k,n)*
72  1 (work(2,k,js)-work(4,k,js))
73  gss(n,m,js) = gss(n,m,js) + p5*tcon(js)*sinnv(k,n)*
74  1 (work(2,k,js)+work(4,k,js))
75  gcc(n,m,js) = gcc(n,m,js) + p5*tcon(js)*cosnv(k,n)*
76  1 (work(1,k,js)+work(3,k,js))
77  END IF
78  END DO
79  END DO
80  END IF
81 !
82 ! INVERSE FOURIER TRANSFORM DE-ALIASED GCON
83 !
84  work(:,:,js) = 0
85 
86  IF(js.GT.1) THEN
87  DO n = 0, ntor
88  DO k = 1, nzeta
89  work(3,k,js) = work(3,k,js) + gcs(n,m,js)*sinnv(k,n)
90  work(4,k,js) = work(4,k,js) + gsc(n,m,js)*cosnv(k,n)
91  IF (.not.lasym) cycle
92  work(1,k,js) = work(1,k,js) + gcc(n,m,js)*cosnv(k,n)
93  work(2,k,js) = work(2,k,js) + gss(n,m,js)*sinnv(k,n)
94  END DO
95  END DO
96  END IF
97 
98  nsmin=tlglob; nsmax=t1rglob
99  DO i = 1, ntheta2
100  DO k = 1, nzeta
101  gcons(k,i,js) = gcons(k,i,js) + (work(3,k,js)*cosmu(i,m)
102  1 + work(4,k,js)*sinmu(i,m))*faccon(m)
103  END DO
104  IF (.not.lasym) cycle
105  DO k = 1, nzeta
106  gcona(k,i,js) = gcona(k,i,js) + (work(1,k,js)*cosmu(i,m)
107  1 + work(2,k,js)*sinmu(i,m))*faccon(m)
108  END DO
109  END DO
110 
111  END DO
112  END DO
113  !END M-LOOP
114 
115  IF (lasym) THEN
116 
117  !EXTEND GCON INTO THETA = PI,2*PI DOMAIN
118  DO js = nsmin, nsmax
119  DO i = 1 + ntheta2, ntheta1
120  ir = ntheta1 + 2 - i
121  DO k = 1, nzeta
122  kk=ireflect_par(k)
123  gcons(k,i,js) = -gcons(kk,ir,js) + gcona(kk,ir,js)
124  END DO
125  END DO
126  END DO
127 
128  !ADD SYMMETRIC, ANTI-SYMMETRIC PIECES IN THETA = 0,PI DOMAIN
129  gcons(:,:ntheta2,nsmin:nsmax) = gcons(:,:ntheta2,nsmin:nsmax)
130  1 + gcona(:,:ntheta2,nsmin:nsmax)
131 
132  END IF
133 
134  DEALLOCATE (work, gcona)
135 
136  CALL second0(taloff)
137  alias_time = alias_time + (taloff-talon)
138 
139  END SUBROUTINE alias_par
140 
141  SUBROUTINE alias(gcons, ztemp, gcs, gsc, gcc, gss)
142  USE vmec_main
143 
144  IMPLICIT NONE
145 C-----------------------------------------------
146 C D u m m y A r g u m e n t s
147 C-----------------------------------------------
148  REAL(dp), PARAMETER :: p5 = 0.5_dp
149  REAL(dp), DIMENSION(ns*nzeta,ntheta3), INTENT(out) :: gcons
150  REAL(dp), DIMENSION(ns*nzeta,ntheta3), INTENT(in) :: ztemp
151  REAL(dp), DIMENSION(ns,0:ntor,0:mpol1) :: gcs, gsc, gcc, gss
152 C-----------------------------------------------
153 C L o c a l V a r i a b l e s
154 C-----------------------------------------------
155  INTEGER :: m, i, ir, jk, jka, n, k, js, l, j
156  REAL(dp), DIMENSION(:,:), ALLOCATABLE :: work, gcona
157 
158 C-----------------------------------------------
159  ALLOCATE (work(ns*nzeta,4), gcona(ns*nzeta,ntheta3))
160 
161  gcons = 0
162  gcona = 0
163 
164  gcs(:,:ntor,:mpol1) = 0; gsc(:,:ntor,:mpol1) = 0
165  gcc(:,:ntor,:mpol1) = 0; gss(:,:ntor,:mpol1) = 0
166 
167  DO m = 1, mpol1 - 1
168  work = 0
169  DO i = 1, ntheta2
170  DO jk = 1, ns*nzeta
171  work(jk,1) = work(jk,1) + ztemp(jk,i)*cosmui(i,m)
172  work(jk,2) = work(jk,2) + ztemp(jk,i)*sinmui(i,m)
173  END DO
174  IF (.not.lasym) cycle
175  !STOP 'CHECK ireflect_par before executing'
176  ir = ntheta1 + 2 - i
177  IF (i .eq. 1) ir = 1
178  DO jk = 1, ns*nzeta
179  jka = ireflect(jk)
180  work(jk,3) = work(jk,3) + ztemp(jka,ir)*cosmui(i,m)
181  work(jk,4) = work(jk,4) + ztemp(jka,ir)*sinmui(i,m)
182  END DO
183  END DO
184 
185  DO n = 0, ntor
186  DO k = 1, nzeta
187  l = ns*(k-1)
188  IF (.not.lasym) THEN
189  DO js = 2,ns
190  gcs(js,n,m) = gcs(js,n,m) + tcon(js)*work(js+l,1)*
191  1 sinnv(k,n)
192  gsc(js,n,m) = gsc(js,n,m) + tcon(js)*work(js+l,2)*
193  1 cosnv(k,n)
194  END DO
195  ELSE
196  DO js = 2,ns
197  gcs(js,n,m) = gcs(js,n,m) + p5*tcon(js)*sinnv(k,n)*
198  1 (work(js+l,1)-work(js+l,3))
199  gsc(js,n,m) = gsc(js,n,m) + p5*tcon(js)*cosnv(k,n)*
200  1 (work(js+l,2)-work(js+l,4))
201  gss(js,n,m) = gss(js,n,m) + p5*tcon(js)*sinnv(k,n)*
202  1 (work(js+l,2)+work(js+l,4))
203  gcc(js,n,m) = gcc(js,n,m) + p5*tcon(js)*cosnv(k,n)*
204  1 (work(js+l,1)+work(js+l,3))
205  END DO
206  END IF
207  END DO
208  END DO
209 !
210 ! INVERSE FOURIER TRANSFORM DE-ALIASED GCON
211 !
212  work = 0
213 
214  DO n = 0, ntor
215  DO k = 1, nzeta
216  l = ns*(k-1)
217  DO js = 2, ns
218  work(js+l,3) = work(js+l,3) + gcs(js,n,m)*sinnv(k,n)
219  work(js+l,4) = work(js+l,4) + gsc(js,n,m)*cosnv(k,n)
220  END DO
221  IF (.not.lasym) cycle
222  DO js = 2, ns
223  work(js+l,1) = work(js+l,1) + gcc(js,n,m)*cosnv(k,n)
224  work(js+l,2) = work(js+l,2) + gss(js,n,m)*sinnv(k,n)
225  END DO
226  END DO
227  END DO
228 
229  DO i = 1, ntheta2
230  DO jk = 1, ns*nzeta
231  gcons(jk,i) = gcons(jk,i) + (work(jk,3)*cosmu(i,m)
232  1 + work(jk,4)*sinmu(i,m))*faccon(m)
233  END DO
234  IF (.not.lasym) cycle
235  DO jk = 1, ns*nzeta
236  gcona(jk,i) = gcona(jk,i) + (work(jk,1)*cosmu(i,m)
237  1 + work(jk,2)*sinmu(i,m))*faccon(m)
238  END DO
239  END DO
240  END DO
241 
242  IF (lasym) THEN
243 
244 ! EXTEND GCON INTO THETA = PI,2*PI DOMAIN
245  DO i = 1 + ntheta2, ntheta1
246  ir = ntheta1 + 2 - i
247  DO jk = 1, ns*nzeta
248  jka = ireflect(jk)
249  gcons(jk,i) = -gcons(jka,ir) + gcona(jka,ir)
250  END DO
251  END DO
252 
253 ! ADD SYMMETRIC, ANTI-SYMMETRIC PIECES IN THETA = 0,PI DOMAIN
254  gcons(:,:ntheta2) = gcons(:,:ntheta2) + gcona(:,:ntheta2)
255 
256  END IF
257 
258  DEALLOCATE (work, gcona)
259 
260  END SUBROUTINE alias