From 7ffdf8e65c90e0048331901588b96715fc319db9 Mon Sep 17 00:00:00 2001 From: JGonzalez Date: Tue, 1 Oct 2024 18:01:27 +0200 Subject: [PATCH] Some updates --- plasmaExpansion.f90 | 21 +++++++++++---------- 1 file changed, 11 insertions(+), 10 deletions(-) diff --git a/plasmaExpansion.f90 b/plasmaExpansion.f90 index b05dc7e..d146046 100644 --- a/plasmaExpansion.f90 +++ b/plasmaExpansion.f90 @@ -349,16 +349,13 @@ program plasmaExpansion L_ref = u_ref * t_ref phi_ref = kb * Temp_ref / qe - ! Set position to calculate cumulative sum of f (non-dimensional units) - rCum = 5.0e-3 / L_ref - ! Set input parameters (remember these have to be in non-dimensional units) Temp0 = 60.0_dp * eV_to_K / Temp_ref - TempF = 60.0_dp * eV_to_K / Temp_ref + TempF = 10.0_dp * eV_to_K / Temp_ref n_ecr = 1.0e19_dp * cm3_to_m3 / n_ref c_s = sqrt(T_to_Z(Temp0) * gam * Temp0) u_bc0 = c_s!sqrt(Temp0) - u_bcF = u_bc0!sqrt(TempF) + u_bcF = sqrt(TempF) n_bc0 = n_ecr / T_to_Z(Temp0) n_bcF = n_bc0!n_ecr*1.0e-1 / T_to_Z(Temp0) @@ -373,6 +370,9 @@ program plasmaExpansion r(i) = dr * float(i-1) + r0 end do + ! Set position to calculate cumulative sum of f (non-dimensional units) + rCum = 5.0e-3 / L_ref + ! Index for cumulative sum rCum_index = minloc(abs(r - rCum), 1) @@ -395,7 +395,7 @@ program plasmaExpansion t0 = 0.0_dp tf = 1.0e-6_dp / t_ref ! tf = 1.0e1_dp * (rf - r0) / c_s - dt = 1.0e-1_dp*dr/(10.0_dp*c_s) + dt = 1.0e-2_dp*dr/(10.0_dp*c_s) nt = nint((tf - t0) / dt) dt = (tf - t0) / float(nt) @@ -500,7 +500,8 @@ program plasmaExpansion n_bc = n_bcF end if call writeOutputBoundary(t, dt, n_bc, u_bc, Temp_bc, T_to_Z(Temp_bc)) - f0(j0:nv) = n_bc / sqrt(PI*Temp_bc) * exp(-(v(j0:nv) - u_bc)**2 / Temp_bc) + f0(j0:nv) = v(j0:nv)**2 / sqrt(PI*Temp_bc**3) * exp(-(v(j0:nv) - u_bc)**2 / Temp_bc) + f0 = f0 * n_bc / (sum(f0)*dv) ! Boundary conditions ! r = r0, v>0 @@ -528,10 +529,10 @@ program plasmaExpansion r(i-1)**2*f_i_old(i-1, j0:nv)) end if - n_i(i) = sum(f_i(i,:)*dv) + n_i(i) = sum(f_i(i,:))*dv if (n_i(i) > 1.0e-10_dp) then - u_i(i) = sum(v(:) *f_i(i,:)*dv) / n_i(i) - E_i(i) = sum(v(:)**2*f_i(i,:)*dv) / n_i(i) + u_i(i) = sum(v(:) *f_i(i,:))*dv / n_i(i) + E_i(i) = sum(v(:)**2*f_i(i,:))*dv / n_i(i) T_i(i) = 2.0_dp*E_i(i) - 2.0_dp*u_i(i)**2 Zave(i) = T_to_Z(T_i(i))