V3FIT
totzsp_mod.f
1  MODULE totzsp_mod
2  USE vmec_main
3  USE timer_sub
4  IMPLICIT NONE
5 
6  INTEGER, PARAMETER, PRIVATE :: m0=0, m1=1, n0=0
7  REAL(dp), ALLOCATABLE, PRIVATE :: work1(:,:,:), work2(:,:)
8  REAL(dp), PRIVATE :: cosmux, sinmux
9 
10  CONTAINS
11 
12 !-------------------------------------------------------------------------------
46 !-------------------------------------------------------------------------------
47  SUBROUTINE totzsps_par(rzl_array, r11, ru1, rv1, z11, zu1, zv1,
48  & lu1, lv1, rcn1, zcn1, ier_flag)
49  USE vmec_params, ONLY: jmin1, jlam, ntmax, rcc, rss, zsc, zcs,
50  & r01_bad_value_flag
51  USE precon2d, ONLY: ictrl_prec2d
52  USE parallel_include_module
53 !-----------------------------------------------
54 ! D u m m y A r g u m e n t s
55 !-----------------------------------------------
56  REAL(dp), DIMENSION(0:ntor,0:mpol1,ns,3*ntmax),
57  & TARGET, INTENT(INOUT) :: rzl_array
58  REAL(dp), DIMENSION(nzeta,ntheta3,ns,0:1),
59  & INTENT(out) :: r11, ru1,
60  & rv1, z11, zu1, zv1, lu1, lv1, rcn1, zcn1
61  INTEGER, INTENT(inout) :: ier_flag
62 !-----------------------------------------------
63 ! L o c a l V a r i a b l e s
64 !-----------------------------------------------
65  INTEGER :: n, m, mparity, k, i, l
66  INTEGER :: ioff, joff, mj, ni, nsz
67  INTEGER :: nsmin, nsmax, js
68  REAL(dp), DIMENSION(:,:,:), POINTER ::
69  & rmncc, rmnss, zmncs, zmnsc, lmncs, lmnsc
70  REAL(dp) :: tbroadon, tbroadoff
71 !-----------------------------------------------
72  CALL second0(tffton)
73 
74  nsmin = t1lglob
75  nsmax = t1rglob
76 
77  rmncc=>rzl_array(:,:,:,rcc) !COS(mu) COS(nv)
78  zmnsc=>rzl_array(:,:,:,zsc+ntmax) !SIN(mu) COS(nv)
79  lmnsc=>rzl_array(:,:,:,zsc+2*ntmax) !SIN(mu) COS(nv)
80  IF (lthreed) THEN
81  rmnss=>rzl_array(:,:,:,rss) !SIN(mu) SIN(nv)
82  zmncs=>rzl_array(:,:,:,zcs+ntmax) !COS(mu) SIN(nv)
83  lmncs=>rzl_array(:,:,:,zcs+2*ntmax) !COS(mu) SIN(nv)
84  END IF
85  rzl_array(:,m1,1,:) = rzl_array(:,m1,2,:)
86 
87  ioff = lbound(rmncc,1)
88  joff = lbound(rmncc,2)
89  IF (lthreed) THEN
90  CALL convert_sym_par(rmnss(:,m1+joff,:), zmncs(:,m1+joff,:),
91  & nsmin, nsmax)
92  END IF
93 
94 !
95 ! ORIGIN EXTRAPOLATION OF M=0 MODES FOR LAMBDA
96 !
97  IF (lthreed .AND. jlam(m0) .GT. 1) THEN
98  lmncs(:,m0+joff,1) = lmncs(:,m0+joff,2)
99  END IF
100 
101 !
102 ! EVOLVE CHIPS BY FORCES IN TOMNSPS WHEN NCURR=1, ICTRL_PREC2D != 0
103 !
104 !SPH071017
105 #if defined(CHI_FORCE)
106  IF (ncurr .EQ. 1) THEN
107  IF (ictrl_prec2d .EQ. 2) THEN
108  lmnsc(n0+ioff,m0+joff,nsmin:nsmax) = chips(nsmin:nsmax)
109  ELSE IF (ictrl_prec2d .NE. 0) THEN
110  chips(nsmin:nsmax) = lmnsc(n0+ioff,m0+joff,nsmin:nsmax)
111  END IF
112  END IF
113 #endif
114 
115  ALLOCATE (work1(nzeta,12,nsmin:nsmax), stat=i)
116  IF (i .ne. 0) THEN
117  stop 'Allocation error in VMEC2000 totzsps'
118  END IF
119 
120  DO js = nsmin, nsmax
121  r11(:,:,js,:) = 0
122  ru1(:,:,js,:) = 0
123  rv1(:,:,js,:) = 0
124  rcn1(:,:,js,:) = 0
125  zcn1(:,:,js,:) = 0
126  z11(:,:,js,:) = 0
127  zu1(:,:,js,:) = 0
128  zv1(:,:,js,:) = 0
129  lu1(:,:,js,:) = 0
130  lv1(:,:,js,:) = 0
131  DO m = 0, mpol1
132  mparity = mod(m,2)
133  mj = m + joff
134  work1(:,:,js) = 0
135 !
136 ! INVERSE TRANSFORM IN N-ZETA, FOR FIXED M
137 !
138  DO n = 0, ntor
139  ni = n + ioff
140  DO k = 1, nzeta
141  work1(k,1,js) = work1(k,1,js)
142  & + rmncc(ni,mj,js)*cosnv(k,n)
143  work1(k,6,js) = work1(k,6,js)
144  & + zmnsc(ni,mj,js)*cosnv(k,n)
145  work1(k,10,js) = work1(k,10,js)
146  & + lmnsc(ni,mj,js)*cosnv(k,n)
147 
148  IF (.NOT.lthreed) cycle
149 
150  work1(k,4,js) = work1(k,4,js)
151  & + rmnss(ni,mj,js)*cosnvn(k,n)
152  work1(k,7,js) = work1(k,7,js)
153  & + zmncs(ni,mj,js)*cosnvn(k,n)
154  work1(k,11,js) = work1(k,11,js)
155  & + lmncs(ni,mj,js)*cosnvn(k,n)
156 
157  work1(k,2,js) = work1(k,2,js)
158  & + rmnss(ni,mj,js)*sinnv(k,n)
159  work1(k,5,js) = work1(k,5,js)
160  & + zmncs(ni,mj,js)*sinnv(k,n)
161  work1(k,9,js) = work1(k,9,js)
162  & + lmncs(ni,mj,js)*sinnv(k,n)
163 
164  work1(k,3,js) = work1(k,3,js)
165  & + rmncc(ni,mj,js)*sinnvn(k,n)
166  work1(k,8,js) = work1(k,8,js)
167  & + zmnsc(ni,mj,js)*sinnvn(k,n)
168  work1(k,12,js) = work1(k,12,js)
169  & + lmnsc(ni,mj,js)*sinnvn(k,n)
170  END DO
171  END DO
172 
173 !
174 ! INVERSE TRANSFORM IN M-THETA, FOR ALL RADIAL, ZETA VALUES
175 !
176  l = 0
177  DO i = 1, ntheta2
178  cosmux = xmpq(m,1)*cosmu(i,m)
179  sinmux = xmpq(m,1)*sinmu(i,m)
180 
181  r11(:,i,js,mparity) = r11(:,i,js,mparity)
182  & + work1(:,1,js)*cosmu(i,m)
183  ru1(:,i,js,mparity) = ru1(:,i,js,mparity)
184  & + work1(:,1,js)*sinmum(i,m)
185  rcn1(:,i,js,mparity) = rcn1(:,i,js,mparity)
186  & + work1(:,1,js)*cosmux
187 
188  z11(:,i,js,mparity) = z11(:,i,js,mparity)
189  & + work1(:,6,js)*sinmu(i,m)
190  zu1(:,i,js,mparity) = zu1(:,i,js,mparity)
191  & + work1(:,6,js)*cosmum(i,m)
192  zcn1(:,i,js,mparity) = zcn1(:,i,js,mparity)
193  & + work1(:,6,js)*sinmux
194 
195  lu1(:,i,js,mparity) = lu1(:,i,js,mparity)
196  & + work1(:,10,js)*cosmum(i,m)
197 
198  IF (.not.lthreed) cycle
199 
200  r11(:,i,js,mparity) = r11(:,i,js,mparity)
201  & + work1(:,2,js)*sinmu(i,m)
202  ru1(:,i,js,mparity) = ru1(:,i,js,mparity)
203  & + work1(:,2,js)*cosmum(i,m)
204  rcn1(:,i,js,mparity) = rcn1(:,i,js,mparity)
205  & + work1(:,2,js)*sinmux
206 
207  rv1(:,i,js,mparity) = rv1(:,i,js,mparity)
208  & + work1(:,3,js)*cosmu(i,m)
209  & + work1(:,4,js)*sinmu(i,m)
210  z11(:,i,js,mparity) = z11(:,i,js,mparity)
211  & + work1(:,5,js)*cosmu(i,m)
212 
213  zu1(:,i,js,mparity) = zu1(:,i,js,mparity)
214  & + work1(:,5,js)*sinmum(i,m)
215  zcn1(:,i,js,mparity) = zcn1(:,i,js,mparity)
216  & + work1(:,5,js)*cosmux
217  zv1(:,i,js,mparity) = zv1(:,i,js,mparity)
218  & + work1(:,7,js)*cosmu(i,m)
219  & + work1(:,8,js)*sinmu(i,m)
220 
221  lu1(:,i,js,mparity) = lu1(:,i,js,mparity)
222  & + work1(:,9,js)*sinmum(i,m)
223  lv1(:,i,js,mparity) = lv1(:,i,js,mparity)
224  & - (work1(:,11,js)*cosmu(i,m)
225  & + work1(:,12,js)*sinmu(i,m))
226  END DO
227  END DO
228  END DO
229 
230  DEALLOCATE (work1)
231 
232  z01(nsmin:nsmax) = zmnsc(n0+ioff,m1+joff,nsmin:nsmax)
233  r01(nsmin:nsmax) = rmncc(n0+ioff,m1+joff,nsmin:nsmax)
234  IF (lactive) THEN
235  IF (rank.EQ.0 .AND. r01(1).EQ.zero) THEN
236  ier_flag = r01_bad_value_flag
237  ELSE IF (rank.EQ.0 .AND. r01(1).NE.zero) THEN
238  dkappa = z01(1)/r01(1)
239  END IF
240  CALL second0(tbroadon)
241  CALL mpi_bcast(dkappa,1, mpi_real8,0,ns_comm,mpi_err)
242  CALL second0(tbroadoff)
243  broadcast_time = broadcast_time + (tbroadoff - tbroadon)
244  END IF
245 
246  CALL second0(tfftoff)
247  totzsps_time = totzsps_time + (tfftoff - tffton)
248  timer(tfft) = timer(tfft) + (tfftoff - tffton)
249 
250  END SUBROUTINE totzsps_par
251 
252 !-------------------------------------------------------------------------------
286 !-------------------------------------------------------------------------------
287  SUBROUTINE totzspa_par(rzl_array, r11, ru1, rv1, z11, zu1, zv1,
288  1 lu1, lv1, rcn1, zcn1)
289  USE vmec_params, ONLY: jmin1, jlam, ntmax, rcs, rsc, zcc, zss
290  USE precon2d, ONLY: ictrl_prec2d
291  USE parallel_include_module
292 C-----------------------------------------------
293 C D u m m y A r g u m e n t s
294 C-----------------------------------------------
295  REAL(dp), DIMENSION(0:ntor,0:mpol1,ns,3*ntmax),
296  1 TARGET, INTENT(inout) :: rzl_array
297  REAL(dp), DIMENSION(nzeta,ntheta3,ns,0:1), INTENT(out) ::
298  1 r11, ru1, rv1, z11, zu1, zv1, lu1, lv1, rcn1, zcn1
299 C-----------------------------------------------
300 C L o c a l V a r i a b l e s
301 C-----------------------------------------------
302  INTEGER :: m, n, mparity, k, i, l, j1
303  INTEGER :: ioff, joff, mj, ni
304  INTEGER :: nsmin, nsmax, js
305  REAL(dp), DIMENSION(:,:,:), POINTER ::
306  1 rmncs, rmnsc, zmncc, zmnss, lmncc, lmnss
307 C-----------------------------------------------
308  CALL second0(tffton)
309  nsmin = t1lglob
310  nsmax = t1rglob
311 
312  rmnsc => rzl_array(:,:,:,rsc) !!SIN(mu) COS(nv)
313  zmncc => rzl_array(:,:,:,zcc+ntmax) !!COS(mu) COS(nv)
314  lmncc => rzl_array(:,:,:,zcc+2*ntmax) !!COS(mu) COS(nv)
315  IF (lthreed) THEN
316  rmncs => rzl_array(:,:,:,rcs) !!COS(mu) SIN(nv)
317  zmnss => rzl_array(:,:,:,zss+ntmax) !!SIN(mu) SIN(nv)
318  lmnss => rzl_array(:,:,:,zss+2*ntmax) !!SIN(mu) SIN(nv)
319  END IF
320 
321 !
322 ! CONVERT FROM INTERNAL XC REPRESENTATION FOR m=1 MODES, R+(at rsc) = .5(rsc + zcc),
323 ! R-(at zcc) = .5(rsc - zcc), TO REQUIRED rsc, zcc FORMS
324 !
325  ioff = lbound(rmnsc,1)
326  joff = lbound(rmnsc,2)
327  CALL convert_asym_par(rmnsc(:,m1+joff,:), zmncc(:,m1+joff,:),
328  & nsmin, nsmax)
329 
330  z00b = zmncc(ioff,joff,ns)
331 
332  ALLOCATE (work1(nzeta,12,nsmin:nsmax), stat=i)
333  IF (i .NE. 0) THEN
334  stop 'Allocation error in VMEC totzspa'
335  END IF
336 
337 !
338 ! INITIALIZATION BLOCK
339 !
340 
341  IF (jlam(m0) .gt. 1) THEN
342  lmncc(:,m0+joff,1) = lmncc(:,m0+joff,2)
343  END IF
344 
345  DO js = nsmin, nsmax
346  r11(:,:,js,:) = 0
347  ru1(:,:,js,:) = 0
348  rv1(:,:,js,:) = 0
349  rcn1(:,:,js,:) = 0
350  zcn1(:,:,js,:) = 0
351  z11(:,:,js,:) = 0
352  zu1(:,:,js,:) = 0
353  zv1(:,:,js,:) = 0
354  lu1(:,:,js,:) = 0
355  lv1(:,:,js,:) = 0
356  DO m = 0, mpol1
357  mparity = mod(m,2)
358  mj = m+joff
359  work1(:,:,js) = 0
360  j1 = jmin1(m)
361 
362  DO n = 0, ntor
363  ni = n+ioff
364  DO k = 1, nzeta
365  work1(k,1,js) = work1(k,1,js)
366  & + rmnsc(ni,mj,js)*cosnv(k,n)
367  work1(k,6,js) = work1(k,6,js)
368  & + zmncc(ni,mj,js)*cosnv(k,n)
369  work1(k,10,js) = work1(k,10,js)
370  & + lmncc(ni,mj,js)*cosnv(k,n)
371 
372  IF (.NOT.lthreed) cycle
373 
374  work1(k,2,js) = work1(k,2,js)
375  & + rmncs(ni,mj,js)*sinnv(k,n)
376  work1(k,3,js) = work1(k,3,js)
377  & + rmnsc(ni,mj,js)*sinnvn(k,n)
378  work1(k,4,js) = work1(k,4,js)
379  & + rmncs(ni,mj,js)*cosnvn(k,n)
380  work1(k,5,js) = work1(k,5,js)
381  & + zmnss(ni,mj,js)*sinnv(k,n)
382  work1(k,7,js) = work1(k,7,js)
383  & + zmnss(ni,mj,js)*cosnvn(k,n)
384  work1(k,8,js) = work1(k,8,js)
385  & + zmncc(ni,mj,js)*sinnvn(k,n)
386  work1(k,9,js) = work1(k,9,js)
387  & + lmnss(ni,mj,js)*sinnv(k,n)
388  work1(k,11,js) = work1(k,11,js)
389  & + lmnss(ni,mj,js)*cosnvn(k,n)
390  work1(k,12,js) = work1(k,12,js)
391  & + lmncc(ni,mj,js)*sinnvn(k,n)
392  END DO
393  END DO
394 
395 !
396 ! INVERSE TRANSFORM IN M-THETA
397 !
398  DO i = 1, ntheta2
399  cosmux = xmpq(m,1)*cosmu(i,m)
400  sinmux = xmpq(m,1)*sinmu(i,m)
401 
402  r11(:,i,js,mparity) = r11(:,i,js,mparity)
403  & + work1(:,1,js)*sinmu(i,m)
404  ru1(:,i,js,mparity) = ru1(:,i,js,mparity)
405  & + work1(:,1,js)*cosmum(i,m)
406  z11(:,i,js,mparity) = z11(:,i,js,mparity)
407  & + work1(:,6,js)*cosmu(i,m)
408  zu1(:,i,js,mparity) = zu1(:,i,js,mparity)
409  & + work1(:,6,js)*sinmum(i,m)
410  lu1(:,i,js,mparity) = lu1(:,i,js,mparity)
411  & + work1(:,10,js)*sinmum(i,m)
412  rcn1(:,i,js,mparity) = rcn1(:,i,js,mparity)
413  & + work1(:,1,js)*sinmux
414  zcn1(:,i,js,mparity) = zcn1(:,i,js,mparity)
415  & + work1(:,6,js)*cosmux
416 
417  IF (.not.lthreed) cycle
418 
419  r11(:,i,js,mparity) = r11(:,i,js,mparity)
420  & + work1(:,2,js)*cosmu(i,m)
421  ru1(:,i,js,mparity) = ru1(:,i,js,mparity)
422  & + work1(:,2,js)*sinmum(i,m)
423  z11(:,i,js,mparity) = z11(:,i,js,mparity)
424  & + work1(:,5,js)*sinmu(i,m)
425  zu1(:,i,js,mparity) = zu1(:,i,js,mparity)
426  & + work1(:,5,js)*cosmum(i,m)
427  lu1(:,i,js,mparity) = lu1(:,i,js,mparity)
428  & + work1(:,9,js)*cosmum(i,m)
429  rcn1(:,i,js,mparity) = rcn1(:,i,js,mparity)
430  & + work1(:,2,js)*cosmux
431  zcn1(:,i,js,mparity) = zcn1(:,i,js,mparity)
432  & + work1(:,5,js)*sinmux
433  rv1(:,i,js,mparity) = rv1(:,i,js,mparity)
434  & + work1(:,3,js)*sinmu(i,m)
435  & + work1(:,4,js)*cosmu(i,m)
436  zv1(:,i,js,mparity) = zv1(:,i,js,mparity)
437  & + work1(:,7,js)*sinmu(i,m)
438  & + work1(:,8,js)*cosmu(i,m)
439  lv1(:,i,js,mparity) = lv1(:,i,js,mparity)
440  & - work1(:,11,js)*sinmu(i,m)
441  & - work1(:,12,js)*cosmu(i,m)
442  END DO
443  END DO
444  END DO
445 
446  DEALLOCATE (work1)
447 
448  CALL second0(tfftoff)
449  totzspa_time = totzspa_time + (tfftoff - tffton)
450  timer(tfft) = timer(tfft) + (tfftoff - tffton)
451 
452  END SUBROUTINE totzspa_par
453 
454  SUBROUTINE convert_sym_par(rmnss, zmncs, nsmin, nsmax)
455 C-----------------------------------------------
456 C D u m m y A r g u m e n t s
457 C-----------------------------------------------
458  INTEGER, INTENT(IN) :: nsmin, nsmax
459  REAL(dp), DIMENSION(0:ntor,ns), INTENT(INOUT) :: rmnss, zmncs
460 C-----------------------------------------------
461 C L o c a l V a r i a b l e s
462 C-----------------------------------------------
463  REAL(dp), DIMENSION(0:ntor,nsmin:nsmax) :: temp
464 C-----------------------------------------------
465 !
466 ! CONVERT FROM INTERNAL REPRESENTATION TO "PHYSICAL" RMNSS, ZMNCS FOURIER FORM
467 ! (for lconm1, rss = zmncs)
468 !
469  IF (lconm1) THEN
470  temp(:,nsmin:nsmax) = rmnss(:,nsmin:nsmax)
471  rmnss(:,nsmin:nsmax) = temp(:,nsmin:nsmax)
472  & + zmncs(:,nsmin:nsmax)
473  zmncs(:,nsmin:nsmax) = temp(:,nsmin:nsmax)
474  & - zmncs(:,nsmin:nsmax)
475  END IF
476 
477  END SUBROUTINE convert_sym_par
478 
479  SUBROUTINE convert_asym_par(rmnsc, zmncc, nsmin, nsmax)
480 C-----------------------------------------------
481 C D u m m y A r g u m e n t s
482 C-----------------------------------------------
483  INTEGER, INTENT(IN) :: nsmin, nsmax
484  REAL(dp), DIMENSION(0:ntor,ns), INTENT(INOUT) :: rmnsc, zmncc
485 C-----------------------------------------------
486 C L o c a l V a r i a b l e s
487 C-----------------------------------------------
488  REAL(dp), DIMENSION(0:ntor,nsmin:nsmax) :: temp
489 C-----------------------------------------------
490 !
491 ! CONVERT FROM INTERNAL REPRESENTATION TO RMNSC, ZMNCC FOURIER FORM
492 !
493  IF (lconm1) THEN
494  temp(:,nsmin:nsmax) = rmnsc(:,nsmin:nsmax)
495  rmnsc(:,nsmin:nsmax) = temp(:,nsmin:nsmax)
496  & + zmncc(:,nsmin:nsmax)
497  zmncc(:,nsmin:nsmax) = temp(:,nsmin:nsmax)
498  & - zmncc(:,nsmin:nsmax)
499  END IF
500 
501  END SUBROUTINE convert_asym_par
502 
503  SUBROUTINE totzsps(rzl_array, r11, ru1, rv1, z11, zu1, zv1,
504  & lu1, lv1, rcn1, zcn1)
505  USE vmec_params, ONLY: jmin1, jlam, ntmax, rcc, rss, zsc, zcs
506  USE precon2d, ONLY: ictrl_prec2d
507 
508 !-----------------------------------------------
509 ! D u m m y A r g u m e n t s
510 !-----------------------------------------------
511  REAL(dp), DIMENSION(ns,0:ntor,0:mpol1,3*ntmax),
512  & TARGET, INTENT(INOUT) :: rzl_array
513  REAL(dp), DIMENSION(ns*nzeta*ntheta3,0:1),
514  & INTENT(OUT) :: r11, ru1,
515  & rv1, z11, zu1, zv1, lu1, lv1, rcn1, zcn1
516 !-----------------------------------------------
517 ! L o c a l V a r i a b l e s
518 !-----------------------------------------------
519  INTEGER :: n, m, mparity, k, i, j1, l, j1l, nsl
520  INTEGER :: ioff, joff, mj, ni, nsz
521  REAL(dp), DIMENSION(:,:,:), POINTER ::
522  & rmncc, rmnss, zmncs, zmnsc, lmncs, lmnsc
523 !-----------------------------------------------
524 !
525 ! WORK1 Array of inverse transforms in toroidal angle (zeta), for all radial positions
526 ! NOTE: ORDERING OF LAST INDEX IS DIFFERENT HERE THAN IN PREVIOUS VMEC2000 VERSIONS
527 !
528 ! CONVERT FROM INTERNAL XC REPRESENTATION FOR m=1 MODES, R+(stored at rss) = .5(rss + zcs),
529 ! R-(stored at zcs) = .5(rss - zcs), TO EXTERNAL ("PHYSICAL") rss, zcs FORMS. NEED THIS EVEN
530 ! WHEN COMPUTING HESSIAN FOR FREE BOUNDARY (rmnss, zmncs at JS=NS needed in vacuum call)
531 !
532  CALL second0(tffton)
533 
534  rmncc => rzl_array(:,:,:,rcc) !!COS(mu) COS(nv)
535  zmnsc => rzl_array(:,:,:,zsc+ntmax) !!SIN(mu) COS(nv)
536  lmnsc => rzl_array(:,:,:,zsc+2*ntmax) !!SIN(mu) COS(nv)
537  IF (lthreed) THEN
538  rmnss => rzl_array(:,:,:,rss) !!SIN(mu) SIN(nv)
539  zmncs => rzl_array(:,:,:,zcs+ntmax) !!COS(mu) SIN(nv)
540  lmncs => rzl_array(:,:,:,zcs+2*ntmax) !!COS(mu) SIN(nv)
541  END IF
542 
543  ioff = lbound(rmncc,2)
544  joff = lbound(rmncc,3)
545 #ifndef _HBANGLE
546  IF (lthreed) THEN
547  CALL convert_sym(rmnss(:,:,m1+joff), zmncs(:,:,m1+joff))
548  END IF
549 #endif
550 
551 !v8.50: Norm for preconditioned R,Z forces: scale to boundary value only
552 !v8.51 Restore hs dependence (1:ns, not just ns)
553 ! fnorm1 = one/SUM(rzl_array(1:ns,:,m1:,1:2*ntmax)**2)
554 
555 !
556 ! ORIGIN EXTRAPOLATION (JS=1) FOR M=1 MODES
557 ! NOTE: PREVIOUS VERSIONS OF VMEC USED TWO-POINT EXTRAPOLATION
558 ! FOR R,Z. HOWEVER,THIS CAN NOT BE USED TO COMPUTE THE
559 ! TRI-DIAG 2D PRECONDITIONER
560 !
561  rzl_array(1,:,m1,:) = rzl_array(2,:,m1,:)
562 
563 !
564 ! ORIGIN EXTRAPOLATION OF M=0 MODES FOR LAMBDA
565 !
566  IF (lthreed .and. jlam(m0) .gt. 1) THEN
567  lmncs(1,:,m0+joff) = lmncs(2,:,m0+joff)
568  END IF
569 
570 !
571 ! EVOLVE CHIPS BY FORCES IN TOMNSPS WHEN NCURR=1, ICTRL_PREC2D != 0
572 !
573 #if defined(CHI_FORCE)
574  IF (ncurr .EQ. 1) THEN
575  IF (ictrl_prec2d .EQ. 2) THEN
576  lmnsc(2:ns,n0+ioff,m0+joff) = chips(2:ns)
577  ELSE IF (ictrl_prec2d .NE. 0) THEN
578  chips(2:ns) = lmnsc(2:ns,n0+ioff,m0+joff)
579  END IF
580  END IF
581 #endif
582  nsz = ns*nzeta
583  ALLOCATE (work2(nsz,12), stat=i)
584  IF (i .NE. 0) THEN
585  stop 'Allocation error in VMEC2000 totzsps'
586  END IF
587 
588  r11 = 0
589  ru1 = 0
590  rv1 = 0
591  rcn1 = 0
592  z11 = 0
593  zu1 = 0
594  zv1 = 0
595  zcn1 = 0
596  lu1 = 0
597  lv1 = 0
598 
599 !
600 ! COMPUTE R, Z, AND LAMBDA IN REAL SPACE
601 ! NOTE: LU = d(Lam)/du, LV = -d(Lam)/dv
602 !
603 
604  DO m = 0, mpol1
605  mparity = mod(m,2)
606  mj = m + joff
607  work2 = 0
608  j1 = jmin1(m)
609 !
610 ! INVERSE TRANSFORM IN N-ZETA, FOR FIXED M
611 !
612  DO n = 0, ntor
613  ni = n+ioff
614  DO k = 1, nzeta
615  l = ns*(k - 1)
616  j1l = j1 + l
617  nsl = ns + l
618  work2(j1l:nsl,1) = work2(j1l:nsl,1)
619  1 + rmncc(j1:ns,ni,mj)*cosnv(k,n)
620  work2(j1l:nsl,6) = work2(j1l:nsl,6)
621  1 + zmnsc(j1:ns,ni,mj)*cosnv(k,n)
622  work2(j1l:nsl,10) = work2(j1l:nsl,10)
623  1 + lmnsc(j1:ns,ni,mj)*cosnv(k,n)
624 
625  IF (.not.lthreed) cycle
626 
627  work2(j1l:nsl,4) = work2(j1l:nsl,4)
628  1 + rmnss(j1:ns,ni,mj)*cosnvn(k,n)
629  work2(j1l:nsl,7) = work2(j1l:nsl,7)
630  1 + zmncs(j1:ns,ni,mj)*cosnvn(k,n)
631  work2(j1l:nsl,11) = work2(j1l:nsl,11)
632  1 + lmncs(j1:ns,ni,mj)*cosnvn(k,n)
633 
634  work2(j1l:nsl,2) = work2(j1l:nsl,2)
635  1 + rmnss(j1:ns,ni,mj)*sinnv(k,n)
636  work2(j1l:nsl,5) = work2(j1l:nsl,5)
637  1 + zmncs(j1:ns,ni,mj)*sinnv(k,n)
638  work2(j1l:nsl,9) = work2(j1l:nsl,9)
639  1 + lmncs(j1:ns,ni,mj)*sinnv(k,n)
640 
641  work2(j1l:nsl,3) = work2(j1l:nsl,3)
642  1 + rmncc(j1:ns,ni,mj)*sinnvn(k,n)
643  work2(j1l:nsl,8) = work2(j1l:nsl,8)
644  1 + zmnsc(j1:ns,ni,mj)*sinnvn(k,n)
645  work2(j1l:nsl,12) = work2(j1l:nsl,12)
646  1 + lmnsc(j1:ns,ni,mj)*sinnvn(k,n)
647  END DO
648  END DO
649 !
650 ! INVERSE TRANSFORM IN M-THETA, FOR ALL RADIAL, ZETA VALUES
651 !
652  l = 0
653  DO i = 1, ntheta2
654  j1l = l + 1
655  nsl = nsz + l
656  l = l + nsz
657  cosmux = xmpq(m,1)*cosmu(i,m)
658  sinmux = xmpq(m,1)*sinmu(i,m)
659 
660  r11(j1l:nsl,mparity) = r11(j1l:nsl,mparity)
661  1 + work2(1:nsz,1)*cosmu(i,m)
662  ru1(j1l:nsl,mparity) = ru1(j1l:nsl,mparity)
663  1 + work2(1:nsz,1)*sinmum(i,m)
664 #ifndef _HBANGLE
665  rcn1(j1l:nsl,mparity) = rcn1(j1l:nsl,mparity)
666  1 + work2(1:nsz,1)*cosmux
667 #endif
668 
669  z11(j1l:nsl,mparity) = z11(j1l:nsl,mparity)
670  1 + work2(1:nsz,6)*sinmu(i,m)
671  zu1(j1l:nsl,mparity) = zu1(j1l:nsl,mparity)
672  1 + work2(1:nsz,6)*cosmum(i,m)
673 #ifndef _HBANGLE
674  zcn1(j1l:nsl,mparity) = zcn1(j1l:nsl,mparity)
675  1 + work2(1:nsz,6)*sinmux
676 #endif
677 
678  lu1(j1l:nsl,mparity) = lu1(j1l:nsl,mparity)
679  1 + work2(1:nsz,10)*cosmum(i,m)
680 
681  IF (.not.lthreed) cycle
682 
683  r11(j1l:nsl,mparity) = r11(j1l:nsl,mparity)
684  1 + work2(1:nsz,2)*sinmu(i,m)
685  ru1(j1l:nsl,mparity) = ru1(j1l:nsl,mparity)
686  1 + work2(1:nsz,2)*cosmum(i,m)
687 #ifndef _HBANGLE
688  rcn1(j1l:nsl,mparity) = rcn1(j1l:nsl,mparity)
689  1 + work2(1:nsz,2)*sinmux
690 #endif
691 
692  rv1(j1l:nsl,mparity) = rv1(j1l:nsl,mparity)
693  1 + work2(1:nsz,3)*cosmu(i,m)
694  1 + work2(1:nsz,4)*sinmu(i,m)
695  z11(j1l:nsl,mparity) = z11(j1l:nsl,mparity)
696  1 + work2(1:nsz,5)*cosmu(i,m)
697 
698  zu1(j1l:nsl,mparity) = zu1(j1l:nsl,mparity)
699  1 + work2(1:nsz,5)*sinmum(i,m)
700 #ifndef _HBANGLE
701  zcn1(j1l:nsl,mparity) = zcn1(j1l:nsl,mparity)
702  1 + work2(1:nsz,5)*cosmux
703 #endif
704  zv1(j1l:nsl,mparity) = zv1(j1l:nsl,mparity)
705  1 + work2(1:nsz,7)*cosmu(i,m)
706  1 + work2(1:nsz,8)*sinmu(i,m)
707 
708  lu1(j1l:nsl,mparity) = lu1(j1l:nsl,mparity)
709  1 + work2(1:nsz,9)*sinmum(i,m)
710  lv1(j1l:nsl,mparity) = lv1(j1l:nsl,mparity)
711  1 - work2(1:nsz,11)*cosmu(i,m)
712  1 - work2(1:nsz,12)*sinmu(i,m)
713  END DO
714  END DO
715 
716  DEALLOCATE (work2)
717 
718  z01(1:ns) = zmnsc(1:ns,n0+ioff,m1+joff)
719  r01(1:ns) = rmncc(1:ns,n0+ioff,m1+joff)
720  IF (r01(1) .eq. zero) THEN
721  stop 'r01(0) = 0 in totzsps_SPH'
722  END IF
723  dkappa = z01(1)/r01(1)
724 
725  END SUBROUTINE totzsps
726 
727  SUBROUTINE convert_sym(rmnss, zmncs)
728 C-----------------------------------------------
729 C D u m m y A r g u m e n t s
730 C-----------------------------------------------
731  REAL(dp), DIMENSION(ns,0:ntor), INTENT(INOUT) :: rmnss, zmncs
732 C-----------------------------------------------
733 C L o c a l V a r i a b l e s
734 C-----------------------------------------------
735  REAL(dp), DIMENSION(ns,0:ntor) :: temp
736 C-----------------------------------------------
737 !
738 ! CONVERT FROM INTERNAL REPRESENTATION TO "PHYSICAL" RMNSS, ZMNCS FOURIER FORM
739 !
740  IF (lconm1) THEN
741  temp = rmnss
742  rmnss = temp + zmncs
743  zmncs = temp - zmncs
744  END IF
745 
746  END SUBROUTINE convert_sym
747 
748 
749  SUBROUTINE totzspa(rzl_array, r11, ru1, rv1, z11, zu1, zv1, lu1,
750  & lv1, rcn1, zcn1)
751  USE vmec_params, ONLY: jmin1, jlam, ntmax, rcs, rsc, zcc, zss
752  USE precon2d, ONLY: ictrl_prec2d
753 
754 C-----------------------------------------------
755 C D u m m y A r g u m e n t s
756 C-----------------------------------------------
757  REAL(dp), DIMENSION(ns,0:ntor,0:mpol1,3*ntmax),
758  1 TARGET, INTENT(inout) :: rzl_array
759  REAL(dp), DIMENSION(ns*nzeta,ntheta3,0:1), INTENT(out) ::
760  1 r11, ru1, rv1, z11, zu1, zv1, lu1, lv1, rcn1, zcn1
761 C-----------------------------------------------
762 C L o c a l V a r i a b l e s
763 C-----------------------------------------------
764  INTEGER :: m, n, mparity, k, i, l, j1, j1l, nsl
765  INTEGER :: ioff, joff, mj, ni
766  REAL(dp), DIMENSION(:,:,:), POINTER ::
767  1 rmncs, rmnsc, zmncc, zmnss, lmncc, lmnss
768 C-----------------------------------------------
769  CALL second0(tffton)
770  rmnsc => rzl_array(:,:,:,rsc) !!SIN(mu) COS(nv)
771  zmncc => rzl_array(:,:,:,zcc+ntmax) !!COS(mu) COS(nv)
772  lmncc => rzl_array(:,:,:,zcc+2*ntmax) !!COS(mu) COS(nv)
773  IF (lthreed) THEN
774  rmncs => rzl_array(:,:,:,rcs) !!COS(mu) SIN(nv)
775  zmnss => rzl_array(:,:,:,zss+ntmax) !!SIN(mu) SIN(nv)
776  lmnss => rzl_array(:,:,:,zss+2*ntmax) !!SIN(mu) SIN(nv)
777  END IF
778 
779 !
780 ! CONVERT FROM INTERNAL XC REPRESENTATION FOR m=1 MODES, R+(at rsc) = .5(rsc + zcc),
781 ! R-(at zcc) = .5(rsc - zcc), TO REQUIRED rsc, zcc FORMS
782 !
783  ioff = lbound(rmnsc,2)
784  joff = lbound(rmnsc,3)
785  CALL convert_asym(rmnsc(:,:,m1+joff), zmncc(:,:,m1+joff))
786 
787  z00b = zmncc(ns,ioff,joff)
788 
789  IF (jlam(m0) .GT. 1) THEN
790  lmncc(1,:,m0+joff) = lmncc(2,:,m0+joff)
791  END IF
792 
793 ! IF (ictrl_prec2d .eq. 3) RETURN
794 
795  ALLOCATE (work2(ns*nzeta,12), stat=i)
796  IF (i .ne. 0) stop 'Allocation error in VMEC totzspa'
797 
798 !
799 ! INITIALIZATION BLOCK
800 !
801  r11 = 0
802  ru1 = 0
803  rv1 = 0
804  z11 = 0
805  zu1 = 0
806  zv1 = 0
807  lu1 = 0
808  lv1 = 0
809  rcn1 = 0
810  zcn1 = 0
811 
812  DO m = 0, mpol1
813  mparity = mod(m,2)
814  mj = m + joff
815  work2 = 0
816  j1 = jmin1(m)
817  DO n = 0, ntor
818  ni = n + ioff
819  DO k = 1, nzeta
820  l = ns*(k - 1)
821  j1l = j1 + l
822  nsl = ns + l
823  work2(j1l:nsl,1) = work2(j1l:nsl,1)
824  1 + rmnsc(j1:ns,ni,mj)*cosnv(k,n)
825  work2(j1l:nsl,6) = work2(j1l:nsl,6)
826  1 + zmncc(j1:ns,ni,mj)*cosnv(k,n)
827  work2(j1l:nsl,10) = work2(j1l:nsl,10)
828  1 + lmncc(j1:ns,ni,mj)*cosnv(k,n)
829 
830  IF (.not.lthreed) cycle
831 
832  work2(j1l:nsl,2) = work2(j1l:nsl,2)
833  1 + rmncs(j1:ns,ni,mj)*sinnv(k,n)
834  work2(j1l:nsl,3) = work2(j1l:nsl,3)
835  1 + rmnsc(j1:ns,ni,mj)*sinnvn(k,n)
836  work2(j1l:nsl,4) = work2(j1l:nsl,4)
837  1 + rmncs(j1:ns,ni,mj)*cosnvn(k,n)
838  work2(j1l:nsl,5) = work2(j1l:nsl,5)
839  1 + zmnss(j1:ns,ni,mj)*sinnv(k,n)
840  work2(j1l:nsl,7) = work2(j1l:nsl,7)
841  1 + zmnss(j1:ns,ni,mj)*cosnvn(k,n)
842  work2(j1l:nsl,8) = work2(j1l:nsl,8)
843  1 + zmncc(j1:ns,ni,mj)*sinnvn(k,n)
844  work2(j1l:nsl,9) = work2(j1l:nsl,9)
845  1 + lmnss(j1:ns,ni,mj)*sinnv(k,n)
846  work2(j1l:nsl,11) = work2(j1l:nsl,11)
847  1 + lmnss(j1:ns,ni,mj)*cosnvn(k,n)
848  work2(j1l:nsl,12) = work2(j1l:nsl,12)
849  1 + lmncc(j1:ns,ni,mj)*sinnvn(k,n)
850  END DO
851  END DO
852 
853 !
854 ! INVERSE TRANSFORM IN M-THETA
855 !
856  DO i = 1, ntheta2
857  cosmux = xmpq(m,1)*cosmu(i,m)
858  sinmux = xmpq(m,1)*sinmu(i,m)
859  r11(:,i,mparity) = r11(:,i,mparity) + work2(:,1)*sinmu(i,m)
860  ru1(:,i,mparity) = ru1(:,i,mparity) + work2(:,1)*cosmum(i,m)
861  z11(:,i,mparity) = z11(:,i,mparity) + work2(:,6)*cosmu(i,m)
862  zu1(:,i,mparity) = zu1(:,i,mparity) + work2(:,6)*sinmum(i,m)
863  lu1(:,i,mparity) = lu1(:,i,mparity)
864  & + work2(:,10)*sinmum(i,m)
865  rcn1(:,i,mparity) = rcn1(:,i,mparity)
866  & + work2(:,1)*sinmux
867  zcn1(:,i,mparity) = zcn1(:,i,mparity)
868  & + work2(:,6)*cosmux
869 
870  IF (.not.lthreed) cycle
871 
872  r11(:,i,mparity) = r11(:,i,mparity)
873  & + work2(:,2)*cosmu(i,m)
874  ru1(:,i,mparity) = ru1(:,i,mparity)
875  & + work2(:,2)*sinmum(i,m)
876  z11(:,i,mparity) = z11(:,i,mparity)
877  & + work2(:,5)*sinmu(i,m)
878  zu1(:,i,mparity) = zu1(:,i,mparity)
879  & + work2(:,5)*cosmum(i,m)
880  lu1(:,i,mparity) = lu1(:,i,mparity)
881  & + work2(:,9)*cosmum(i,m)
882  rcn1(:,i,mparity) = rcn1(:,i,mparity)
883  & + work2(:,2)*cosmux
884  zcn1(:,i,mparity) = zcn1(:,i,mparity)
885  & + work2(:,5)*sinmux
886  rv1(:,i,mparity) = rv1(:,i,mparity)
887  & + work2(:,3)*sinmu(i,m)
888  & + work2(:,4)*cosmu(i,m)
889  zv1(:,i,mparity) = zv1(:,i,mparity)
890  & + work2(:,7)*sinmu(i,m)
891  & + work2(:,8)*cosmu(i,m)
892  lv1(:,i,mparity) = lv1(:,i,mparity)
893  & - work2(:,11)*sinmu(i,m)
894  & - work2(:,12)*cosmu(i,m)
895  END DO
896  END DO
897 
898  DEALLOCATE (work2)
899 
900  END SUBROUTINE totzspa
901 
902 
903  SUBROUTINE convert_asym(rmnsc, zmncc)
904 C-----------------------------------------------
905 C D u m m y A r g u m e n t s
906 C-----------------------------------------------
907  REAL(dp), DIMENSION(ns,0:ntor), INTENT(INOUT) :: rmnsc, zmncc
908 C-----------------------------------------------
909 C L o c a l V a r i a b l e s
910 C-----------------------------------------------
911  REAL(dp), DIMENSION(ns,0:ntor) :: temp
912 C-----------------------------------------------
913 !
914 ! CONVERT FROM INTERNAL REPRESENTATION TO RMNSC, ZMNCC FOURIER FORM
915 !
916  IF (lconm1) THEN
917  temp = rmnsc
918  rmnsc = temp + zmncc
919  zmncc = temp - zmncc
920  END IF
921 
922  END SUBROUTINE convert_asym
923 
924  END MODULE totzsp_mod
Generated on Thu Mar 5 2020 15:49:26 for V3FIT by   doxygen 1.8.17