First attempt at self-consistent field solve. There are many glitches on... (e8f03111) · Commits · Kumar, Atul / LinearFokkerPlanck Axisymmetric

InputFiles/xp_TestDrag.in

+6 −6

Original line number	Diff line number	Diff line
		&params_nml
		! Simulation name:
		! ===============
		params%fileDescriptor = '2021_03_28b',
		params%fileDescriptor = '2021_03_28d',

		! Magnetic field input data:
		! =========================
		@@ -12,7 +12,7 @@ params%nz = 501,

		! Simulation conditions:
		! =====================
		params%NC = 200000, ! Total number of particles
		params%NC = 500000, ! Total number of particles
		params%NS = 20000, ! Total number of time steps
		params%dt = 0.5E-7, ! Time step in [s]
		params%G = 2.5E+20, ! Real particle injection rate
		@@ -43,8 +43,8 @@ params%iDrag = .true.,

		! Collision operator conditions:
		! ==============================
		params%Te0 = 250., ! eV
		params%Ti0 = 250., ! eV
		params%Te0 = 600., ! eV
		params%Ti0 = 600., ! eV
		params%ne0 = 0.5E+19, ! m**-3
		params%Aion = 1., ! Ion mass
		params%Zion = 1., ! Charge number for main ion species
		@@ -54,11 +54,11 @@ params%CollOperType = 2,

		! Particle Boundary Conditions:
		! ============================
		params%BC_Type = 2, ! 1: Isotropic plasma source, 2: NBI, 3: Periodic
		params%BC_Type = 1, ! 1: Isotropic plasma source, 2: NBI, 3: Periodic
		params%BC_zp_mean = 0.0, ! Mean particle injection
		params%BC_zp_std = 0.3, ! STD of particle injection
		params%BC_Ep = 2000., ! Drift kinetic energy
		params%BC_Tp = 10, ! Thermal kinetic energy
		params%BC_Tp = 600, ! Thermal kinetic energy
		params%BC_xip = 0.707, ! Mean pitch angle where xip = cos(theta) = vpar/v

		! Particle Initial Conditions:

src/Fields.f90

+55 −5

Original line number	Diff line number	Diff line
		SUBROUTINE AdvanceEfield(plasma,mesh,params)
		SUBROUTINE AdvanceEfield(mesh,params)
		USE LOCAL
		USE dataTYP
		USE PhysicalConstants
		@@ -6,7 +6,6 @@ USE PhysicalConstants
		IMPLICIT NONE

		! Declare interface variables:
		TYPE(plasmaTYP), INTENT(IN) :: plasma
		TYPE(meshTYP), INTENT(INOUT) :: mesh
		TYPE(paramsTYP), INTENT(IN) :: params

		@@ -25,7 +24,7 @@ dz = mesh%dzm
		Te = params%Te0

		! Plasma density gradient:
		rng = (\ (i, i = 3, (NZ + 2)) \)
		rng = (/ (i, i = 3, (NZ + 2)) /)
		dn = (n(rng+1) - n(rng-1))/(2*dz)

		! Electric field:
		@@ -34,10 +33,61 @@ E(rng) = -Te*dn/n(rng)

		! Apply smoothing:
		frame = 9
		CALL MovingMean(E(rng),NZ,frame)
		CALL MovingMean_F(E,NZg,frame)

		! Output:
		mesh%E = E

		RETURN
		END SUBROUTINE AdvanceEfield

		! =======================================================================================================
		SUBROUTINE MovingMean_F(y,NX,k)
		! =======================================================================================================
		USE LOCAL
		USE dataTYP

		IMPLICIT NONE

		! Define interface variables:
		REAL(r8) , DIMENSION(NX), INTENT(INOUT) :: y
		INTEGER(i4), INTENT(IN) :: NX
		INTEGER(i4), INTENT(INOUT) :: k

		! Define local variables:
		REAL(r8), DIMENSION(NX) :: ym
		REAL(r8) :: yd
		INTEGER(i4) :: s, ii, jj, istart, iend, N

		! Check frame:
		IF (MOD(k,2) .EQ. 0) THEN
		k = k + 1
		END IF

		! Half frame size:
		s = (k-1)/2

		DO ii = 1,NX
		! Start and end of frame:
		istart = ii-s
		iend = ii+s

		! Correct frame at edges:
		IF (istart .LE. 0 ) istart = 1
		IF (iend .GT. NX) iend = NX

		! Effective frame size:
		N = iend - istart + 1

		! Calculate mean:
		yd = 0.
		DO jj = istart,iend
		yd = yd + y(jj)/N
		END DO
		ym(ii) = yd
		END DO

		y = ym

		END SUBROUTINE MovingMean_F

src/Modules.f90

+8 −0

Original line number	Diff line number	Diff line
		@@ -790,6 +790,14 @@ SUBROUTINE SaveData(output,dir1)
		OPEN(unit=8,file=fileName,form="unformatted",status="unknown")
		WRITE(8) output%U
		CLOSE(unit=8)
		fileName = trim(trim(dir1)//'/'//'ddB_mesh.out')
		OPEN(unit=8,file=fileName,form="unformatted",status="unknown")
		WRITE(8) output%ddB
		CLOSE(unit=8)
		fileName = trim(trim(dir1)//'/'//'E_mesh.out')
		OPEN(unit=8,file=fileName,form="unformatted",status="unknown")
		WRITE(8) output%E
		CLOSE(unit=8)


		! Saving pcount to file:

src/MoveParticlePack.f90

+7 −1

Original line number	Diff line number	Diff line
		@@ -487,14 +487,20 @@ USE dataTYP
		IMPLICIT NONE

		! Define interface variables:
		INTEGER(i4), INTENT(IN) :: NX, k
		REAL(r8) , DIMENSION(NX), INTENT(INOUT) :: y
		INTEGER(i4), INTENT(IN) :: NX
		INTEGER(i4), INTENT(INOUT) :: k

		! Define local variables:
		REAL(r8), DIMENSION(NX) :: ym
		REAL(r8) :: yd
		INTEGER(i4) :: s, ii, jj, istart, iend, N

		! Check frame:
		IF (MOD(k,2) .EQ. 0) THEN
		k = k + 1
		END IF

		! Half frame size:
		s = (k-1)/2

src/linFP.f90

+3 −2

Original line number	Diff line number	Diff line
		@@ -183,7 +183,7 @@ NS_loop: DO j = 1,params%NS
		!$OMP END PARALLEL DO

		! Field solve:
		CALL AdvanceEfield(plasma,mesh,params)
		CALL AdvanceEfield(mesh,params)

		! Interpolate electromagnetic fields to particle positions:
		CALL InterpolateElectromagneticFields(plasma,mesh,params)
		@@ -248,7 +248,8 @@ NS_loop: DO j = 1,params%NS
		output%Tpar(:,k) = mesh%Tpar
		output%Tper(:,k) = mesh%Tper
		output%U(:,k) = mesh%U

		output%ddB(:,k) = mesh%ddB
		output%E(:,k) = mesh%E
		! Increment counter:
		k = k + 1
		END IF