4 INTEGER,
PARAMETER :: nvertex = 4
9 INTEGER :: xpt, ypt, phi_plane
13 REAL(rprec) :: br, bp, bz
18 TYPE (MAGFIELD),
POINTER :: pbfld
24 TYPE (VERTEX),
POINTER :: pprev
25 TYPE (PBFIELD) :: bvector(nvertex)
26 LOGICAL :: lalloc(nvertex)
30 TYPE (VERTEX),
POINTER :: pverx
34 REAL(rprec) :: hR, hZ, hP
35 INTEGER,
ALLOCATABLE :: ibs_cnt(:)
36 TYPE (PVERTEX),
ALLOCATABLE,
TARGET :: tail(:)
37 TYPE (VERTEX),
POINTER :: root
38 TYPE (VERTEX),
POINTER :: new_root, prev_vertex
39 LOGICAL :: lsearch(nvertex)
43 SUBROUTINE initialize_list (number_lists, hR_in, hZ_in, hP_in)
45 INTEGER,
INTENT(in) :: number_lists
46 REAL(rprec),
INTENT(in) :: hR_in, hZ_in, hP_in
50 IF (nlist .le. 0) stop
' Must specify nlist > 0'
52 ALLOCATE (tail(nlist), ibs_cnt(nlist))
55 NULLIFY (tail(ilist)%pverx)
60 hr = hr_in; hz = hz_in; hp = hp_in
62 END SUBROUTINE initialize_list
65 SUBROUTINE add_vertex (rin, zin, phi_pt, list_index, bfield)
67 TYPE (VERTEX),
POINTER :: new_vertex, swap
68 TYPE (MAGFIELD),
POINTER :: MagFld
69 INTEGER,
INTENT(in) :: phi_pt, list_index
70 REAL(rprec),
INTENT(in) :: rin, zin
71 REAL(rprec) :: Rreal, Zreal, PhiReal
72 INTEGER :: rpt, zpt, xpt, ypt
81 IF (list_index .gt. nlist) stop
'list_index > nlist in add_vertex'
82 IF (list_index .le. 0) stop
'list_index <= 0 in add_vertex'
88 root => tail(list_index)%pverx
91 IF (
ASSOCIATED(root))
THEN
92 new_vertex%pprev => root
94 NULLIFY(new_vertex%pprev)
97 new_vertex%p1 =
point(rpt, zpt, phi_pt)
100 CALL find_vertex (new_vertex, lexists)
106 DEALLOCATE (new_vertex)
107 IF (
ASSOCIATED(prev_vertex))
THEN
117 DO icount = 1, nvertex
118 IF (.not.lsearch(icount))
THEN
120 ALLOCATE (new_vertex%bvector(icount)%pbfld)
121 new_vertex%lalloc(icount) = .true.
122 magfld => new_vertex%bvector(icount)%pbfld
125 IF (icount.eq.2 .or. icount.eq.4) xpt = xpt+1
126 IF (icount.ge.3) ypt = ypt+1
127 rreal = xpt*hr; zreal = ypt*hz
128 phireal = (phi_pt-1)*hp
129 CALL bfield (rreal, phireal, zreal,
130 1 magfld%br, magfld%bp, magfld%bz)
131 ibs_cnt(list_index) = ibs_cnt(list_index) + 1
133 new_vertex%lalloc(icount) = .false.
138 tail(list_index)%pverx => root
140 END SUBROUTINE add_vertex
143 SUBROUTINE find_vertex (new_node, lexists)
145 TYPE (VERTEX),
POINTER :: new_node
146 LOGICAL,
INTENT(out) :: lexists
147 TYPE (VERTEX),
POINTER :: stepit, previous
148 INTEGER :: xpt(nvertex), ypt(nvertex), xvert,
149 1 yvert, xstep, ystep, icount, ivertex
156 xpt(1) = new_node%p1%xpt; ypt(1) = new_node%p1%ypt
157 xpt(2) = xpt(1)+1; ypt(2) = ypt(1)
158 xpt(3) = xpt(1); ypt(3) = ypt(1)+1
159 xpt(4) = xpt(2); ypt(4) = ypt(3)
162 NULLIFY(previous);
NULLIFY(new_root);
NULLIFY(prev_vertex)
165 lsearch(1:nvertex) = .false.
168 DO WHILE (icount.lt.nvertex .and.
ASSOCIATED(stepit))
169 xstep = stepit%p1%xpt; ystep = stepit%p1%ypt
174 IF (xstep.eq.xpt(1) .and. ystep.eq.ypt(1))
THEN
176 IF (
ASSOCIATED(previous)) previous%pprev => stepit%pprev
177 prev_vertex => previous
181 DO ivertex = 1, nvertex
182 xvert = xpt(ivertex); yvert = ypt(ivertex)
184 IF ((xstep.gt.xvert) .or. (ystep.gt.yvert) .or.
185 1 (xstep.lt.xvert-1) .or. (ystep.lt.yvert-1) .or.
186 2 lsearch(ivertex) ) cycle
188 SELECT CASE (ivertex)
190 IF (xstep.eq.xvert)
THEN
191 new_node%bvector(1) = stepit%bvector(3)
193 IF (ystep .eq. yvert)
THEN
194 new_node%bvector(1) = stepit%bvector(2)
196 new_node%bvector(1) = stepit%bvector(4)
201 IF (xstep.lt.xvert)
THEN
202 new_node%bvector(2) = stepit%bvector(4)
204 IF (ystep .eq. yvert)
THEN
205 new_node%bvector(2) = stepit%bvector(1)
207 new_node%bvector(2) = stepit%bvector(3)
212 IF (xstep .eq. xvert)
THEN
213 new_node%bvector(3) = stepit%bvector(1)
215 IF (ystep .eq. yvert)
THEN
216 new_node%bvector(3) = stepit%bvector(2)
218 new_node%bvector(3) = stepit%bvector(4)
223 IF (xstep .lt. xvert)
THEN
224 new_node%bvector(4) = stepit%bvector(2)
226 IF (ystep .eq. yvert)
THEN
227 new_node%bvector(4) = stepit%bvector(1)
229 new_node%bvector(4) = stepit%bvector(3)
236 IF (icount .ge. nvertex)
EXIT
237 lsearch(ivertex) = .true.
242 stepit => stepit%pprev
246 END SUBROUTINE find_vertex
249 SUBROUTINE b_at_point (rin, zin, br, bp, bz)
254 REAL(rprec),
INTENT(in) :: rin, zin
255 REAL(rprec),
INTENT(out) :: br, bp, bz
256 REAL(rprec),
DIMENSION(nvertex) :: brvac, bpvac, bzvac
260 TYPE (MAGFIELD),
POINTER :: bvac
261 REAL(rprec),
PARAMETER :: zero = 0, one = 1
263 REAL(rprec) :: ri, zj, dri, dzj, fint(nvertex)
294 bvac => root%bvector(i)%pbfld
304 fint(3) = dzj - fint(4)
305 fint(2) = dri - fint(4)
306 fint(1) = 1 + fint(4) - (dri + dzj)
308 IF (dri.lt.zero .or. dri.gt.one .or.
309 1 dzj.lt.zero .or. dzj.gt.one)
THEN
310 WRITE (6,
'(2(a,1pe12.4))')
311 1
'Possible interp. error in b_from_coil: dri = ',
319 END SUBROUTINE b_at_point
321 SUBROUTINE cleanup_list
323 TYPE (VERTEX),
POINTER :: stepit
324 INTEGER :: icount, ilist
327 root => tail(ilist)%pverx
328 DO WHILE (
ASSOCIATED(root))
330 DO icount = 1, nvertex
331 IF (root%lalloc(icount))
332 1
DEALLOCATE (root%bvector(icount)%pbfld)
339 DEALLOCATE (tail, ibs_cnt)
341 END SUBROUTINE cleanup_list
343 END MODULE vertex_list