Merge pull request #71 from ganjinzero1/develop

Removing cell crossing near the X-point for a HFS 1D flux tube

bbb/geometry.m

@@ -790,7 +790,7 @@ call set_indirect_address(0)
 
 c...  For 1D poloidal flux-tube geometry (nysol=1), reset cell vertex
 c...  locations near ix=ixpt2 on iy=1 surface if nxpt2msh or nxpt2psh > 1
-      if (nxpt2msh+nxpt2psh > 0) call reset1dmeshpt
+      if (nxpt2msh+nxpt2psh > 0 .or. nxpt1msh+nxpt1psh > 0) call reset1dmeshpt
 
 *     Define guard cells around the edge of the mesh --
       call guardc
@@ -2893,10 +2893,14 @@ real xn(0:nx+1,0:ny+1), xvn(0:nx+1,0:ny+1),
       Use(Share)          # nxomit,nxpt2msh,nxpt2psh,rxpt2msh,zxpt2psh
 
 c...  local scalars
-      integer ix,ix2
+      integer ix,ix2,ix1 
       real lenginc,theta
 
-      ix2 = nxomit+ixpt2(1)    		#consider outer flux-surface only
+      if (nxomit .ge. 0) then
+        ix2 = nxomit+ixpt2(1)    		#consider outer flux-surface only
+      else
+	ix2 = ixpt2(1)
+      endif
 
 c...  First fix iy=1 cell vertices above X-point (ixpt2)
       theta = rxpt2msh/zxpt2msh   #defines curvature of flux surface
@@ -2916,6 +2920,32 @@ real xn(0:nx+1,0:ny+1), xvn(0:nx+1,0:ny+1),
         zm(ix+1,1,3) = zm(ix,1,4)
       enddo
 
+#..   Do the inner X-pt
+#..
+#..
+      if (nxomit .ge. 0) then
+        ix1 = nxomit+ixpt1(1)+1
+      else
+        ix1 = ixpt1(1)+1
+      endif
+c...  First fix iy=1 cell vertices above X-point (ixpt2)
+      theta = rxpt1msh/zxpt1msh   #defines curvature of flux surface
+      do ix = ix1, ix1+nxpt1msh-1
+        rm(ix,1,4) = rm(ix,1,3)
+        zm(ix,1,4) = zm(ix,1,3) + zxpt1msh
+        rm(ix-1,1,3) = rm(ix,1,4)
+        zm(ix-1,1,3) = zm(ix,1,4)
+      enddo
+
+c...  Second, fix iy=1 cell vertices below X-point
+      theta = rxpt1psh/zxpt1psh   #defines curvature of flux surface
+      do ix = ix1-1, ix1-nxpt1psh, -1
+        rm(ix,1,3) = rm(ix,1,4)
+        zm(ix,1,3) = zm(ix,1,4) - zxpt1psh
+        rm(ix+1,1,4) = rm(ix,1,3)
+        zm(ix+1,1,4) = zm(ix,1,3)
+      enddo
+
       return
       end
 c***  End of subroutine reset1dmeshpt ****

com/com.v

@@ -276,6 +276,12 @@ zxpt2msh        real /0.01/ [m] # spacing iy=1 1D vertices reset above ixpt2
 rxpt2msh        real /0.0/  [m] # spacing iy=1 1D vertices reset above ixpt2
 zxpt2psh        real /0.01/ [m] # spacing iy=1 1D vertices reset below ixpt2
 rxpt2psh        real /0.0/  [m] # spacing iy=1 1D vertices reset below ixpt2
+nxpt1msh        integer /0/     # number of 1D cells above ixpt1 reset/no-cross
+nxpt1psh        integer /0/     # number of 1D cells below ixpt1 reset/no-cross
+zxpt1msh        real /0.01/ [m] # spacing iy=1 1D vertices reset above ixpt1
+rxpt1msh        real /0.0/  [m] # spacing iy=1 1D vertices reset above ixpt1
+zxpt1psh        real /0.01/ [m] # spacing iy=1 1D vertices reset below ixpt1
+rxpt1psh        real /0.0/  [m] # spacing iy=1 1D vertices reset below ixpt1
 ismpsym         integer /0/     +gridgen # =1 re-constructs "guard" cells at midplane
                                 # of "dnbot" via up/down symmetry
 isudsym         integer /0/     +input #=1 up-down symmetric setup (only down part is modeled)