4 INTEGER,
PARAMETER :: ncur = 100, msmax = 64, nresmax = 10
5 INTEGER,
PARAMETER :: iout = 10, jout = 11, kout = 12, lout = 13
6 INTEGER :: nfp, ngroups, ntmax, npmax, lopt, nres, numres
7 INTEGER,
DIMENSION(0:nresmax) :: mres, nimax, nrp
8 REAL(rprec) :: drsurf, residsum, raxiscon
9 REAL(rprec),
DIMENSION(msmax) :: ri, zi, di, drho, iota
10 REAL(rprec),
DIMENSION(0:nresmax) :: vres, rres, dres
11 REAL(rprec),
DIMENSION(0:nresmax) :: trace, deter, resid
13 REAL(rprec),
DIMENSION(0:nresmax) :: rres_min, rres_max, reswt
14 REAL(rprec),
DIMENSION(1:nresmax) :: rmax_isl, rmin_isl
15 REAL(rprec),
DIMENSION(1:nresmax) :: zmax_isl, zmin_isl
16 REAL(rprec),
DIMENSION(1:nresmax) :: shear_r, delu_r
17 REAL(rprec),
DIMENSION(1:nresmax) :: shear_p, delu_p
18 REAL(rprec),
DIMENSION(1:nresmax) :: area_isl, flux_isl
19 REAL(rprec),
DIMENSION(1:nresmax) :: alpha
20 REAL(rprec),
DIMENSION(1:nresmax) :: sigma_resid
21 LOGICAL,
DIMENSION(1:nresmax) :: lsymm
22 REAL(rprec),
DIMENSION(:,:),
ALLOCATABLE :: ropnts, zopnts
23 REAL(rprec),
DIMENSION(:,:),
ALLOCATABLE :: r_isl, z_isl
24 LOGICAL,
DIMENSION(:,:),
ALLOCATABLE :: lli, lri
25 REAL(rprec),
DIMENSION(:,:),
ALLOCATABLE :: rt, zt
26 REAL(rprec),
DIMENSION(msmax) :: rmag, zmag
27 INTEGER :: ierr, mc, msurf, nvar, nstell_coils
28 INTEGER :: ierr_mpi, myid, numprocs
29 REAL(rprec),
DIMENSION(ncur) :: vac_extcur, xc
30 REAL(rprec) :: raxis, zaxis, r0, z0, dsmin, dsurf, tsurf, v0
31 REAL(rprec),
DIMENSION(nresmax) :: atol, rtol
32 LOGICAL :: laxis, lprt, lisize, lsurf, ldisplace, lscreen
33 LOGICAL,
DIMENSION(:),
ALLOCATABLE :: lopnt
36 REAL(rprec) :: as_array(3)
40 TYPE (angles_shifts),
DIMENSION(20) :: angles, angles_lo,
42 TYPE (angles_shifts),
DIMENSION(20) :: shifts, shifts_lo,
43 1 shifts_hi, sigma_shifts_lo, sigma_shifts_hi
44 TYPE (angles_shifts),
DIMENSION(20,10,10) ::
45 1 control_pt_init, control_pt_final
46 CHARACTER(LEN=100) :: extension, coils_file_extension
48 namelist /nlst_vacopt/ vac_extcur, msurf, ntmax, ri, zi, lopt,
49 1 r0, z0, v0, dsmin, dsurf, tsurf, nfp, laxis, raxis, zaxis,
50 2 vres, rres, rres_min, rres_max, reswt, nrp, dres, mres, numres,
51 3 raxiscon, rmin_isl, rmax_isl, zmin_isl, zmax_isl, shear_r,
52 4 shear_p, alpha, lsymm, lisize, lsurf, lprt, ldisplace,
53 5 nstell_coils, nvar, rtol, atol, coils_file_extension,
54 6 rmag, zmag, sigma_resid,
55 7 angles, angles_lo, angles_hi, sigma_angles_lo, sigma_angles_hi,
56 8 shifts, shifts_lo, shifts_hi, sigma_shifts_lo, sigma_shifts_hi,
65 SUBROUTINE initialize_vacfield (arg1)
69 INTEGER :: n, np, istat, invac
70 CHARACTER*(*),
INTENT(in) :: arg1
71 CHARACTER :: invac_file*100
74 extension = trim(arg1)
75 n = index(extension,
".", back=.true.)
76 extension = extension(n+1:)
77 coils_file_extension = trim(extension)
84 CALL safe_open (invac, istat, invac_file,
'old',
'formatted')
85 IF (istat .ne. 0)
THEN
86 print *,
'Error opening namelist input file', invac_file
90 CALL read_vacopt_namelist(invac, istat)
92 IF (istat .ne. 0)
THEN
93 print *,
' Error reading ', trim(invac_file),
95 print *,
' Be sure EXTCUR -> VAC_EXTUR!'
105 CALL initialize_biotsavart (vac_extcur, coils_file_extension,
108 END SUBROUTINE initialize_vacfield
110 SUBROUTINE read_vacopt_namelist (iunit, istat)
111 INTEGER :: iunit, istat
142 DO n = 1,
SIZE(angles,1)
145 sigma_angles_lo(n) =
angles_shifts( (/1._dp, 1._dp, 1._dp/) )
146 sigma_shifts_lo(n) =
angles_shifts( (/1._dp, 1._dp, 1._dp/) )
147 sigma_angles_hi(n) =
angles_shifts( (/1._dp, 1._dp, 1._dp/) )
148 sigma_shifts_hi(n) =
angles_shifts( (/1._dp, 1._dp, 1._dp/) )
152 control_pt_init(n,:,:) =
171 rres_min(0) = 0.70_dp
172 rres_max(0) = 0.73_dp
173 rres_min(1) = 1.33_dp
174 rres_max(1) = 1.36_dp
175 rres_min(2) = 0.88_dp
176 rres_max(2) = 0.94_dp
177 rres_min(3) = 1.37_dp
178 rres_max(3) = 1.40_dp
188 rmin_isl(1) = 1.33_dp
189 rmax_isl(1) = 1.36_dp
190 zmin_isl(1) = -0.40_dp
191 zmax_isl(1) = 0.40_dp
192 rmin_isl(2) = 0.88_dp
193 rmax_isl(2) = 0.94_dp
194 zmin_isl(2) = -0.10_dp
195 zmax_isl(2) = 0.10_dp
196 rmin_isl(3) = 1.37_dp
197 rmax_isl(3) = 1.40_dp
198 zmin_isl(3) = -0.40_dp
199 zmax_isl(3) = 0.40_dp
203 shear_p(1) = -0.3875_dp
204 shear_p(2) = -0.3875_dp
205 shear_p(3) = -0.2813_dp
219 di(n) = dsmin + n*dsurf
220 ri(n) = raxis + di(n)*cos(tsurf)
221 zi(n) = zaxis + di(n)*sin(tsurf)
226 READ (iunit, nml=nlst_vacopt, iostat=istat)
228 END SUBROUTINE read_vacopt_namelist
232 USE write_array_generic
237 INTEGER :: iunit, istat
242 CHARACTER(LEN=100) :: form1, form2
245 WRITE(iunit,
'(a12)')
'&NLST_VACOPT'
246 WRITE(iunit,100)
'LDISPLACE = ', ldisplace
247 WRITE(iunit,100)
'LAXIS = ', laxis
248 WRITE(iunit,100)
'LPRT = ', lprt
249 WRITE(iunit,100)
'LISIZE = ', lisize
250 WRITE(iunit,200)
'LOPT = ', lopt
251 WRITE(iunit,200)
'MSURF = ', msurf
252 WRITE(iunit,200)
'NTMAX = ', ntmax
253 WRITE(iunit,200)
'NFP = ', nfp
254 WRITE(iunit,200)
'NSTELL_COILS = ', nstell_coils
255 WRITE(iunit,200)
'NVAR = ', nvar
256 WRITE(iunit,300)
'RAXIS = ', raxis
257 WRITE(iunit,300)
'R0 = ', r0
258 WRITE(iunit,300)
'Z0 = ', z0
259 WRITE(iunit,300)
'V0 = ', v0
260 WRITE(iunit,300)
'DSURF = ', dsurf
261 WRITE(iunit,300)
'DSMIN = ', dsmin
262 WRITE(iunit,300)
'TSURF = ', tsurf
263 WRITE(iunit,200)
'NUMRES = ', numres
264 WRITE(iunit,
'(1x,3a)')
"COILS_FILE_EXTENSION = '",
265 1 trim(coils_file_extension),
"'"
266 IF (msmax .gt. 0)
THEN
272 IF (nresmax .gt. 0)
THEN
273 CALL write_array(iunit,
'MRES', mres, nresmax+1, low_index=0)
274 CALL write_array(iunit,
'VRES', vres, nresmax+1, low_index=0)
275 CALL write_array(iunit,
'RRES', rres, nresmax+1, low_index=0)
276 CALL write_array(iunit,
'NRP', nrp, nresmax+1, low_index=0)
277 CALL write_array(iunit,
'RRES_MIN', rres_min, nresmax+1,
279 CALL write_array(iunit,
'RRES_MAX', rres_max, nresmax+1,
281 CALL write_array(iunit,
'RESWT', reswt, nresmax+1, low_index=0)
282 CALL write_array(iunit,
'DRES', dres, nresmax+1, low_index=0)
284 IF (ncur .gt. 0)
THEN
285 CALL write_array(iunit,
'VAC_EXTCUR', vac_extcur, ncur)
287 IF (nresmax .gt. 0)
THEN
288 CALL write_array(iunit,
'RMIN_ISL', rmin_isl, nresmax)
289 CALL write_array(iunit,
'RMAX_ISL', rmax_isl, nresmax)
290 CALL write_array(iunit,
'ZMIN_ISL', zmin_isl, nresmax)
291 CALL write_array(iunit,
'ZMAX_ISL', zmax_isl, nresmax)
292 CALL write_array(iunit,
'SHEAR_R', shear_r, nresmax)
293 CALL write_array(iunit,
'SHEAR_P', shear_p, nresmax)
296 CALL write_array(iunit,
'SIGMA_RESID', sigma_resid, nresmax)
299 DO i = 1, nstell_coils
301 form1 =
'(2x,a,i1,a,1p3e15.6)'
302 form2 =
'(2x,a,i1,2(a1,i2.2),a,1p3e15.6)'
304 form1 =
'(2x,a,i2,a,1p3e15.6)'
305 form2 =
'(2x,a,i2,2(a1,i2.2),a,1p3e15.6)'
307 WRITE(iunit, form1, iostat=istat)
308 1
'ANGLES(',i,
') = ', angles(i)
309 WRITE(iunit, form1, iostat=istat)
310 1
'SHIFTS(',i,
') = ', shifts(i)
311 WRITE(iunit, form1, iostat=istat)
312 1
'ANGLES_LO(',i,
') = ', angles_lo(i)
313 WRITE(iunit, form1, iostat=istat)
314 1
'SHIFTS_LO(',i,
') = ', shifts_lo(i)
315 WRITE(iunit, form1, iostat=istat)
316 1
'ANGLES_HI(',i,
') = ', angles_hi(i)
317 WRITE(iunit, form1, iostat=istat)
318 1
'SHIFTS_HI(',i,
') = ', shifts_hi(i)
319 WRITE(iunit, form1, iostat=istat)
320 1
'SIGMA_ANGLES_LO(',i,
') = ', sigma_angles_lo(i)
321 WRITE(iunit, form1, iostat=istat)
322 1
'SIGMA_SHIFTS_LO(',i,
') = ', sigma_shifts_lo(i)
323 WRITE(iunit, form1, iostat=istat)
324 1
'SIGMA_ANGLES_HI(',i,
') = ', sigma_angles_hi(i)
325 WRITE(iunit, form1, iostat=istat)
326 1
'SIGMA_SHIFTS_HI(',i,
') = ', sigma_shifts_hi(i)
328 DO j = 1,
SIZE(control_pt_init,2)
329 DO k = 1,
SIZE(control_pt_init,3)
330 WRITE(iunit, form2, iostat=istat)
331 1
'CONTROL_PT_INIT(',i,
',',j,
',',k,
') = ',
332 2 control_pt_init(i,j,k)
338 WRITE(iunit,
'(a)')
'/'
340 100
FORMAT(4(1x,a,l2,
','))
341 200
FORMAT(4(1x,a,i6,
','))
342 300
FORMAT(3(1x,a,1pe12.4,
','))
343 350
FORMAT(2(1x,a,1pe12.4,
','))
344 450
FORMAT(1x,a,1pe12.4,
',')
348 SUBROUTINE write_invac2 (filename, istat)
350 USE write_array_generic
356 INTEGER :: iunit, istat
357 CHARACTER(LEN=*) :: filename
359 CALL safe_open(iunit, istat, filename,
'replace',
'formatted')
360 IF (istat .ne. 0)
RETURN
366 END SUBROUTINE write_invac2
368 END MODULE vacfield_mod