First step for including the ExB backg. shear flow

Makefile

@@ -118,6 +118,7 @@ FOBJ=$(OBJDIR)/advance_field_mod.o $(OBJDIR)/array_mod.o $(OBJDIR)/auxval.o \
 $(OBJDIR)/basic_mod.o $(OBJDIR)/coeff_mod.o $(OBJDIR)/closure_mod.o $(OBJDIR)/circular_mod.o \
 $(OBJDIR)/collision_mod.o $(OBJDIR)/nonlinear_mod.o $(OBJDIR)/control.o \
 $(OBJDIR)/diagnose.o $(OBJDIR)/diagnostics_par_mod.o $(OBJDIR)/endrun.o \
+$(OBJDIR)/ExB_shear_flow_mod.o \
 $(OBJDIR)/fields_mod.o $(OBJDIR)/fourier_mod.o $(OBJDIR)/geometry_mod.o \
 $(OBJDIR)/ghosts_mod.o $(OBJDIR)/grid_mod.o $(OBJDIR)/inital.o \
 $(OBJDIR)/initial_par_mod.o $(OBJDIR)/lag_interp_mod.o $(OBJDIR)/main.o \
@@ -206,6 +207,11 @@ $(OBJDIR)/time_integration_mod.o $(OBJDIR)/utility_mod.o $(OBJDIR)/CLA_mod.o
  	 $(OBJDIR)/prec_const_mod.o $(OBJDIR)/basic_mod.o
 	$(F90) -c $(F90FLAGS) $(FPPFLAGS) $(EXTMOD) $(EXTINC) src/endrun.F90 -o $@
 
+ $(OBJDIR)/ExB_shear_flow_mod.o : src/ExB_shear_flow_mod.F90 \
+	 $(OBJDIR)/basic_mod.o $(OBJDIR)/grid_mod.o $(OBJDIR)/prec_const_mod.o\
+	 $(OBJDIR)/geometry_mod.o $(OBJDIR)/numerics_mod.o
+	$(F90) -c $(F90FLAGS) $(FPPFLAGS) $(EXTMOD) $(EXTINC) src/ExB_shear_flow_mod.F90 -o $@
+
  $(OBJDIR)/fields_mod.o : src/fields_mod.F90 \
    $(OBJDIR)/prec_const_mod.o
 	$(F90) -c $(F90FLAGS) $(FPPFLAGS) $(EXTMOD) $(EXTINC) src/fields_mod.F90 -o $@

src/ExB_shear_flow_mod.F90

+MODULE ExB_shear_flow
+    ! This module contains the necessary tools to implement ExB shearing flow effects.
+    ! The algorithm is taken from the presentation of Hammett et al. 2006 (APS) and
+    ! it the one used in GS2.
+    USE prec_const, ONLY: xp
+
+    IMPLICIT NONE
+    ! Variables
+    REAL(xp), PUBLIC, PROTECTED :: gamma_E = 0._xp ! ExB background shearing rate \gamma_E
+    REAL(xp), DIMENSION(:), ALLOCATABLE, PUBLIC, PROTECTED :: sky_ExB      ! shift of the kx modes, kx* = kx + s(ky)
+    INTEGER,  DIMENSION(:), ALLOCATABLE, PUBLIC, PROTECTED :: jump_ExB     ! jump to do to shift the kx grids
+    LOGICAL,  DIMENSION(:), ALLOCATABLE, PUBLIC, PROTECTED :: shiftnow_ExB ! Indicates if there is a line to shift
+
+    ! Routines
+    PUBLIC :: Setup_ExB_shear_flow, Update_ExB_shear_flow
+
+CONTAINS
+
+    ! Setup the variables for the ExB shear
+    SUBROUTINE Setup_ExB_shear_flow(g_E)
+        USE basic,      ONLY : time
+        USE grid,       ONLY : local_nky, kyarray, Lx
+        USE prec_const, ONLY : PI
+        IMPLICIT NONE
+        REAL(xp), INTENT(IN) :: g_E ! Input shearing rate (comes from model module)
+        ! local var
+        INTEGER :: iky
+        REAL(xp):: inv_dkx
+
+        ! Store the shearing rate in this module
+        gamma_E = g_E
+
+        ! Setup the ExB shift
+        ALLOCATE(sky_ExB(local_nky))
+        sky_ExB = 0._xp
+
+        ! Setup the jump array and shifting flag
+        ALLOCATE(shiftnow_ExB(local_nky))
+        shiftnow_ExB = .FALSE.
+        ALLOCATE(jump_ExB(local_nky))
+        jump_ExB     = 0
+    END SUBROUTINE Setup_ExB_shear_flow
+
+    ! Update according to the current ExB shear value
+    ! -the grids
+    ! -the spatial operators 
+    ! -the fields by imposing a shift on kx
+    ! Then update the ExB shear value for the next time step
+    SUBROUTINE Update_ExB_shear_flow
+        USE basic,      ONLY : dt
+        USE grid,       ONLY : local_nky, kyarray, update_grids
+        USE prec_const, ONLY : PI
+        USE geometry,   ONLY : gxx,gxy,gyy,inv_hatB2, evaluate_magn_curv
+        USE numerics,   ONLY : evaluate_EM_op, evaluate_kernels
+        IMPLICIT NONE
+        ! local var
+        INTEGER :: iky
+        REAL(xp):: inv_dkx
+
+        ! update the grids
+        CALL update_grids(sky_ExB(iky),gxx,gxy,gyy,inv_hatB2)
+
+        ! update the EM operators and the kernels
+        CALL evaluate_kernels
+        CALL evaluate_EM_op
+
+        ! shift all fields and correct the shift value
+        CALL Shift_fields
+
+        ! update the ExB shift
+        shiftnow_ExB = .FALSE.
+        DO iky = 1,local_Nky
+            sky_ExB(iky)      = sky_ExB(iky) - kyarray(iky)*gamma_E*dt
+            jump_ExB(iky)     = NINT(sky_ExB(iky)*inv_dkx)
+            shiftnow_ExB(iky) = (abs(jump_ExB(iky)) .GT. 0)
+        ENDDO
+        CONTAINS
+
+        SUBROUTINE Shift_fields
+            USE grid,  ONLY: local_nky, kxarray, kyarray, local_nkx, deltakx, &
+                             kx_min, kx_max
+            USE fields,ONLY: moments, phi, psi
+            IMPLICIT NONE
+            ! local var
+            INTEGER :: iky, ikx, ikx_s
+            REAL(xp):: inv_dkx
+            DO iky = 1,local_Nky
+                IF(shiftnow_ExB(iky)) THEN
+                    ! shift all fields
+                    DO ikx = 1,local_nkx
+                        ikx_s = ikx + jump_ExB(iky)
+                        IF( (kxarray(iky,ikx) .GE. kx_min) .AND. (kxarray(iky,ikx) .LE. kx_max) ) THEN
+                            moments(:,:,:,iky,ikx,:,:) = moments(:,:,:,iky,ikx_s,:,:)
+                            phi(iky,ikx,:)             = phi(iky,ikx_s,:)
+                            psi(iky,ikx,:)             = psi(iky,ikx_s,:)
+                        ELSE
+                            moments(:,:,:,iky,ikx,:,:) = 0._xp
+                            phi(iky,ikx,:)             = 0._xp
+                            psi(iky,ikx,:)             = 0._xp
+                        ENDIF
+                    ENDDO
+                    ! correct the shift value s(ky) for this row 
+                    sky_ExB(iky) = sky_ExB(iky) - jump_ExB(iky)*deltakx
+                ENDIF
+            ENDDO
+        END SUBROUTINE Shift_fields
+    END SUBROUTINE Update_ExB_shear_flow
+
+
+END MODULE ExB_shear_flow
\ No newline at end of file