VMEC 8.52
3D Equilibrium solver with nested flux surfaces.
Loading...
Searching...
No Matches
tomnsp.f90
Go to the documentation of this file.
1
3
17SUBROUTINE tomnsps(frzl_array, &
18 armn, brmn, crmn, &
19 azmn, bzmn, czmn, &
20 blmn, clmn, &
21 arcon, azcon )
22 USE vmec_main
23 USE vmec_params, ONLY: jlam, jmin2, ntmax, rcc, rss, zsc, zcs
24
25 use dbgout
26
27 IMPLICIT NONE
28
29 REAL(rprec), DIMENSION(ns,0:ntor,0:mpol1,3*ntmax), TARGET, INTENT(out) :: frzl_array
30 REAL(rprec), DIMENSION(ns*nzeta*ntheta3,0:1), INTENT(in) :: &
31 armn, brmn, crmn, azmn, bzmn, czmn, blmn, clmn, arcon, azcon
32
33 INTEGER :: jmax, m, mparity, i, n, k, l, nsz
34 INTEGER :: ioff, joff, mj, ni, nsl, j2, j2l, jl, jll, jmaxl, js
35 REAL(rprec), DIMENSION(:,:,:), POINTER :: frcc, frss, fzcs, fzsc, flcs, flsc
36 REAL(rprec), ALLOCATABLE, DIMENSION(:,:) :: work1
37 REAL(rprec), DIMENSION(:), ALLOCATABLE :: tempr, tempz
38
39 frcc => frzl_array(:,:,:,rcc) !!COS(mu) COS(nv)
40 fzsc => frzl_array(:,:,:,zsc+ntmax) !!SIN(mu) COS(nv)
41 flsc => frzl_array(:,:,:,zsc+2*ntmax) !!SIN(mu) COS(nv)
42 IF (lthreed) THEN
43 frss => frzl_array(:,:,:,rss) !!SIN(mu) SIN(nv)
44 fzcs => frzl_array(:,:,:,zcs+ntmax) !!COS(mu) SIN(nv)
45 flcs => frzl_array(:,:,:,zcs+2*ntmax) !!COS(mu) SIN(nv)
46 END IF
47
48 nsz = ns*nzeta
49
50 ALLOCATE (work1(nsz,12), tempr(nsz), tempz(nsz), stat=i)
51 IF (i .ne. 0) stop 'Allocation error in VMEC2000 tomnsps'
52
53 ioff = lbound(frcc,2)
54 joff = lbound(frcc,3)
55
56 frzl_array = 0.0_dp
57
58 ! exclude forces on boundary if free-boundary module is not active
59 jmax = ns
60 IF (ivac .lt. 1) jmax = ns1
61
62 ! BEGIN FOURIER TRANSFORM
63 !
64 ! FRmn = ARmn - d(BRmn)/du + d(CRmn)/dv
65 ! FZmn = AZmn - d(BZmn)/du + d(CZmn)/dv
66 ! FLmn = - d(BLmn)/du + d(CLmn)/dv
67 !
68 ! NOTE: sinmumi = -m sin(mu)
69 ! sinnvn = -n*nfp sin(nv)
70 DO m = 0, mpol1
71
72 mparity = mod(m,2)
73 mj = m+joff
74
75 j2 = jmin2(m)
76 jl = jlam(m)
77
78 work1 = 0.0_dp
79
80 ! DO THETA (U) INTEGRATION FIRST ON HALF INTERVAL (0 < U < PI)
81 l = 0
82 DO i = 1, ntheta2
83 jll = l+1 ! start of poloidal slice
84 nsl = l+nsz ! end of poloidal slice
85 l = l+nsz ! jump to next poloidal slice
86
87 tempr(:) = armn(jll:nsl,mparity) &
88! #ifndef _HBANGLE
89 + xmpq(m,1)*arcon(jll:nsl,mparity)
90! #end /* ndef _HBANGLE */
91
92 tempz(:) = azmn(jll:nsl,mparity) &
93! #ifndef _HBANGLE
94 + xmpq(m,1)*azcon(jll:nsl,mparity)
95! #end /* ndef _HBANGLE */
96
97 work1(:, 1) = work1(:, 1) + tempr(:)*cosmui(i,m) + brmn(jll:nsl,mparity)*sinmumi(i,m)
98 work1(:, 7) = work1(:, 7) + tempz(:)*sinmui(i,m) + bzmn(jll:nsl,mparity)*cosmumi(i,m)
99 work1(:,11) = work1(:,11) + blmn(jll:nsl,mparity)*cosmumi(i,m)
100
101 IF (.not.lthreed) cycle
102
103 work1(:, 2) = work1(:, 2) - crmn(jll:nsl,mparity)*cosmui(i,m)
104 work1(:, 3) = work1(:, 3) + tempr(:)*sinmui(i,m) + brmn(jll:nsl,mparity)*cosmumi(i,m)
105 work1(:, 4) = work1(:, 4) - crmn(jll:nsl,mparity)*sinmui(i,m)
106 work1(:, 5) = work1(:, 5) + tempz(:)*cosmui(i,m) + bzmn(jll:nsl,mparity)*sinmumi(i,m)
107 work1(:, 6) = work1(:, 6) - czmn(jll:nsl,mparity)*cosmui(i,m)
108 work1(:, 8) = work1(:, 8) - czmn(jll:nsl,mparity)*sinmui(i,m)
109
110 work1(:, 9) = work1(:, 9) + blmn(jll:nsl,mparity)*sinmumi(i,m)
111 work1(:,10) = work1(:,10) - clmn(jll:nsl,mparity)*cosmui(i,m)
112 work1(:,12) = work1(:,12) - clmn(jll:nsl,mparity)*sinmui(i,m)
113 END DO
114
115 ! NEXT, DO ZETA (V) TRANSFORM
116 DO n = 0, ntor
117 ni = n+ioff
118 l = 0
119 DO k = 1, nzeta
120 j2l = l+j2 ! start of radial slice for R,Z
121 jmaxl = l+jmax ! end of radial slice for R,Z
122
123 jll = l+jl ! start of radial slice for lambda
124 nsl = l+ns ! end of radial slice for lambda
125
126 l = l+ns ! jump to next radial slice
127
128 frcc(j2:jmax,ni,mj) = frcc(j2:jmax,ni,mj) + work1(j2l:jmaxl,1)*cosnv(k,n)
129 fzsc(j2:jmax,ni,mj) = fzsc(j2:jmax,ni,mj) + work1(j2l:jmaxl,7)*cosnv(k,n)
130 flsc(jl:ns,ni,mj) = flsc(jl:ns,ni,mj) + work1(jll:nsl,11)*cosnv(k,n)
131
132 IF (.not.lthreed) cycle
133
134 frcc(j2:jmax,ni,mj) = frcc(j2:jmax,ni,mj) + work1(j2l:jmaxl,2)*sinnvn(k,n)
135 fzsc(j2:jmax,ni,mj) = fzsc(j2:jmax,ni,mj) + work1(j2l:jmaxl,8)*sinnvn(k,n)
136 frss(j2:jmax,ni,mj) = frss(j2:jmax,ni,mj) + work1(j2l:jmaxl,3)*sinnv(k,n) &
137 + work1(j2l:jmaxl,4)*cosnvn(k,n)
138 fzcs(j2:jmax,ni,mj) = fzcs(j2:jmax,ni,mj) + work1(j2l:jmaxl,5)*sinnv(k,n) &
139 + work1(j2l:jmaxl,6)*cosnvn(k,n)
140
141 flsc(jl:ns,ni,mj) = flsc(jl:ns,ni,mj) + work1(jll:nsl,12)*sinnvn(k,n)
142 flcs(jl:ns,ni,mj) = flcs(jl:ns,ni,mj) + work1(jll:nsl, 9)*sinnv(k,n) &
143 + work1(jll:nsl,10)*cosnvn(k,n)
144 END DO
145 END DO
146 END DO
147
148 DEALLOCATE (work1, tempr, tempz)
149
150 if (open_dbg_context("tomnsps", num_eqsolve_retries)) then
151
152 call add_real_3d("frcc", ns, ntor1, mpol, frcc)
153 call add_real_3d("fzsc", ns, ntor1, mpol, fzsc)
154 call add_real_3d("flsc", ns, ntor1, mpol, flsc)
155
156 if (lthreed) then
157 call add_real_3d("frss", ns, ntor1, mpol, frss)
158 call add_real_3d("fzcs", ns, ntor1, mpol, fzcs)
159 call add_real_3d("flcs", ns, ntor1, mpol, flcs)
160 else
161 call add_null("frss")
162 call add_null("fzcs")
163 call add_null("flcs")
164 end if
165
166 call close_dbg_out()
167 end if
168
169END SUBROUTINE tomnsps
170
184SUBROUTINE tomnspa(frzl_array, &
185 armn, brmn, crmn, &
186 azmn, bzmn, czmn, &
187 blmn, clmn, &
188 arcon, azcon )
189 USE vmec_main
190 USE vmec_params, ONLY: jlam, jmin2, ntmax, rsc, rcs, zcc, zss
191
192 use dbgout
193
194 IMPLICIT NONE
195
196 REAL(rprec), DIMENSION(ns,0:ntor,0:mpol1,3*ntmax), TARGET, INTENT(inout) :: frzl_array
197 REAL(rprec), DIMENSION(ns*nzeta,ntheta3,0:1), INTENT(in) :: &
198 armn, brmn, crmn, azmn, bzmn, czmn, blmn, clmn, arcon, azcon
199
200 INTEGER :: jmax, m, mparity, i, n, k, l, nsz
201 INTEGER :: ioff, joff, mj, ni, nsl, j2, j2l, jl, jll, jmaxl, js
202 REAL(rprec), DIMENSION(:,:,:), POINTER :: frcs, frsc, fzcc, fzss, flcc, flss
203 REAL(rprec), DIMENSION(:), ALLOCATABLE :: temp1, temp3
204 REAL(rprec), DIMENSION(:,:), ALLOCATABLE :: work1
205
206 frsc => frzl_array(:,:,:,rsc) !!R-SIN(mu) COS(nv)
207 fzcc => frzl_array(:,:,:,zcc+ntmax) !!Z-COS(mu) COS(nv)
208 flcc => frzl_array(:,:,:,zcc+2*ntmax) !!L-COS(mu) COS(nv)
209 IF (lthreed) THEN
210 frcs => frzl_array(:,:,:,rcs) !!R-COS(mu) SIN(nv)
211 fzss => frzl_array(:,:,:,zss+ntmax) !!Z-SIN(mu) SIN(nv)
212 flss => frzl_array(:,:,:,zss+2*ntmax) !!L-SIN(mu) SIN(nv)
213 END IF
214
215 nsz = ns*nzeta
216 ALLOCATE (work1(nsz,12), temp1(nsz), temp3(nsz), stat=i)
217 IF (i .ne. 0) stop 'Allocation error in VMEC tomnspa'
218
219 ioff = lbound(frsc,2)
220 joff = lbound(frsc,3)
221
222 jmax = ns
223 IF (ivac .lt. 1) jmax = ns1
224
225 ! BEGIN FOURIER TRANSFORM
226 DO m = 0, mpol1
227 mparity = mod(m,2)
228
229 mj = m+joff
230
231 j2 = jmin2(m)
232 jl = jlam(m)
233
234 work1 = 0.0_dp
235
236 ! DO THETA (U) TRANSFORM FIRST
237 DO i = 1, ntheta2
238 temp1(:) = armn(:,i,mparity) &
239! #ifndef _HBANGLE
240 + xmpq(m,1)*arcon(:,i,mparity)
241! #end /* ndef _HBANGLE */
242
243 temp3(:) = azmn(:,i,mparity) &
244! #ifndef _HBANGLE
245 + xmpq(m,1)*azcon(:,i,mparity)
246! #end /* ndef _HBANGLE */
247
248 work1(:,3) = work1(:,3) + temp1(:)*sinmui(i,m) + brmn(:,i,mparity)*cosmumi(i,m)
249 work1(:,5) = work1(:,5) + temp3(:)*cosmui(i,m) + bzmn(:,i,mparity)*sinmumi(i,m)
250 work1(:,9) = work1(:,9) + blmn(:,i,mparity)*sinmumi(i,m)
251
252 IF (.not.lthreed) cycle
253
254 work1(:,1) = work1(:,1) + temp1(:)*cosmui(i,m) + brmn(:,i,mparity)*sinmumi(i,m)
255 work1(:,2) = work1(:,2) - crmn(:,i,mparity)*cosmui(i,m)
256 work1(:,4) = work1(:,4) - crmn(:,i,mparity)*sinmui(i,m)
257 work1(:,6) = work1(:,6) - czmn(:,i,mparity)*cosmui(i,m)
258 work1(:,7) = work1(:,7) + temp3(:)*sinmui(i,m) + bzmn(:,i,mparity)*cosmumi(i,m)
259 work1(:,8) = work1(:,8) - czmn(:,i,mparity)*sinmui(i,m)
260 work1(:,10) = work1(:,10) - clmn(:,i,mparity)*cosmui(i,m)
261 work1(:,11) = work1(:,11) + blmn(:,i,mparity)*cosmumi(i,m)
262 work1(:,12) = work1(:,12) - clmn(:,i,mparity)*sinmui(i,m)
263 END DO
264
265 ! NEXT, DO ZETA (V) TRANSFORM
266 DO n = 0, ntor
267 ni = n+ioff
268 DO k = 1, nzeta
269 l = ns*(k-1) ! current slice offset
270
271 j2l = j2+l ! start of radial slice for R,Z
272 jmaxl = jmax+l ! end of radial slice for R,Z
273
274 jll = jl+l ! start of radial slice for lambda
275 nsl = ns+l ! end of radial slice for lambda
276
277 frsc(j2:jmax,ni,mj) = frsc(j2:jmax,ni,mj) + work1(j2l:jmaxl,3)*cosnv(k,n)
278 fzcc(j2:jmax,ni,mj) = fzcc(j2:jmax,ni,mj) + work1(j2l:jmaxl,5)*cosnv(k,n)
279 flcc(jl:ns,ni,mj) = flcc(jl:ns,ni,mj) + work1(jll:nsl,9)*cosnv(k,n)
280
281 IF (.not.lthreed) cycle
282
283 frsc(j2:jmax,ni,mj) = frsc(j2:jmax,ni,mj) + work1(j2l:jmaxl,4)*sinnvn(k,n)
284 fzcc(j2:jmax,ni,mj) = fzcc(j2:jmax,ni,mj) + work1(j2l:jmaxl,6)*sinnvn(k,n)
285 frcs(j2:jmax,ni,mj) = frcs(j2:jmax,ni,mj) + work1(j2l:jmaxl,1)*sinnv(k,n) &
286 + work1(j2l:jmaxl,2)*cosnvn(k,n)
287 fzss(j2:jmax,ni,mj) = fzss(j2:jmax,ni,mj) + work1(j2l:jmaxl,7)*sinnv(k,n) &
288 + work1(j2l:jmaxl,8)*cosnvn(k,n)
289 flcc(jl:ns,ni,mj) = flcc(jl:ns,ni,mj) + work1(jll:nsl,10)*sinnvn(k,n)
290 flss(jl:ns,ni,mj) = flss(jl:ns,ni,mj) + work1(jll:nsl,11)*sinnv(k,n) &
291 + work1(jll:nsl,12)*cosnvn(k,n)
292 END DO
293 END DO
294 END DO
295
296 DEALLOCATE (work1, temp1, temp3)
297
298 if (open_dbg_context("tomnspa", num_eqsolve_retries)) then
299
300 call add_real_3d("frsc", ns, ntor1, mpol, frsc)
301 call add_real_3d("fzcc", ns, ntor1, mpol, fzcc)
302 call add_real_3d("flcc", ns, ntor1, mpol, flcc)
303
304 if (lthreed) then
305 call add_real_3d("frcs", ns, ntor1, mpol, frcs)
306 call add_real_3d("fzss", ns, ntor1, mpol, fzss)
307 call add_real_3d("flss", ns, ntor1, mpol, flss)
308 else
309 call add_null("frcs")
310 call add_null("fzss")
311 call add_null("flss")
312 end if
313
314 call close_dbg_out()
315 end if
316
317END SUBROUTINE tomnspa
logical function open_dbg_context(context_name, repetition, id)
check if any output is desired for the current iteration check if the given context should be openend...
Definition dbgout.f90:17
integer ivac
counts number of free-boundary iterations
integer num_eqsolve_retries
real(rprec), dimension(0:mpol1d, 3) xmpq
spectral condensation weighting factors
Definition vmec_main.f90:81
logical lthreed
integer, dimension(0:mpold), parameter jmin2
starting js(m) values for which R,Z are evolved
integer, dimension(0:mpold), parameter jlam
starting js(m) values for which Lambda is evolved
integer ntmax
number of contributing Fourier basis function (can be 1, 2 or 4); assigned in read_indata()
subroutine tomnspa(frzl_array, armn, brmn, crmn, azmn, bzmn, czmn, blmn, clmn, arcon, azcon)
Fourier-transform anti-symmetric forces from real space to Fourier space.
Definition tomnsp.f90:189
subroutine tomnsps(frzl_array, armn, brmn, crmn, azmn, bzmn, czmn, blmn, clmn, arcon, azcon)
Fourier-transform symmetric forces from real space to Fourier space.
Definition tomnsp.f90:22