From 0c27b98e2ebb41fb863f3ac0f88b45157fd4ae61 Mon Sep 17 00:00:00 2001 From: Jorge Gonzalez Date: Thu, 10 Apr 2025 08:49:01 +0200 Subject: [PATCH] n_e as boundary condition Now n_e is given as the density at the boundary and n_i at the boundary is calculated once Z is known. This aims to eliminate the iterative process. --- moduleOutput.f90 | 10 +++++----- scripts_python/plotBC.py | 12 +++++++----- scripts_python/plotCumF.py | 6 ++++-- scripts_python/plotMom.py | 11 ++++++----- scripts_python/readBC.py | 9 +++++---- vlaplex.f90 | 28 +++++++++++----------------- 6 files changed, 38 insertions(+), 38 deletions(-) diff --git a/moduleOutput.f90 b/moduleOutput.f90 index 0640145..9a3a882 100644 --- a/moduleOutput.f90 +++ b/moduleOutput.f90 @@ -175,13 +175,13 @@ module output end do end subroutine writeOutputMom - subroutine writeOutputBoundary(t, dt, n, u, Temp, Zinj) + subroutine writeOutputBoundary(t, dt, n, u, Temp, Z_Tne, Zinj) use constantParameters, only: eV_to_K use referenceValues, only: t_ref, n_ref, u_ref, Temp_ref integer, intent(in):: t real(dp), intent(in):: dt - real(dp), intent(in):: n, u, Temp, Zinj + real(dp), intent(in):: n, u, Temp, Z_Tne, Zinj character(len=6), parameter:: filename = 'bc.csv' logical:: res @@ -189,14 +189,14 @@ module output if (.not. res) then write (*, '(A, A)') 'Writing: ', filename open(unit=dataBC_id, file=pathOutput // filename, action='write', position='append') - write(dataBC_id, '(A,4(' // formatSep // ',A))') 't (s)', 'n (m^-3)', 'u (m s^-1)', 'T (eV)','Zinj' + write(dataBC_id, '(A,5(' // formatSep // ',A))') 't (s)', 'n_i (m^-3)', 'u_i (m s^-1)', 'T_i (eV)', 'Zinj','Z_Tne' close(dataBC_id) end if open(unit=dataBC_id, file=pathOutput // filename, action='write', position='append') - write(dataBC_id, '(' // formatFloat // ',4('// formatSep // ',' // formatFloat // '))') & - t*dt*t_ref, n*n_ref, u*u_ref, Temp*Temp_ref/eV_to_K, Zinj + write(dataBC_id, '(' // formatFloat // ',5('// formatSep // ',' // formatFloat // '))') & + t*dt*t_ref, n*n_ref, u*u_ref, Temp*Temp_ref/eV_to_K, Zinj, Z_Tne close(dataBC_id) diff --git a/scripts_python/plotBC.py b/scripts_python/plotBC.py index 445f3b9..207869e 100644 --- a/scripts_python/plotBC.py +++ b/scripts_python/plotBC.py @@ -4,15 +4,17 @@ import numpy as np import readBC -paths = ['../polytropic_80ns_T30/'] -time, n, u, T, Zinj = readBC.read(paths[0] + 'bc.csv') +paths = ['../2025-04-10_08.40.09/'] +time, n_i, u_i, T_i, Zinj, Z_Tne = readBC.read(paths[0] + 'bc.csv') fig, ax = plt.subplots() +n_e = n_i * Zinj -plt.plot(time, n / n[0], label = f"$n_i$ ($\\times {n[0] * 1e-6:.0e} \\; cm^{{-3}})$") -plt.plot(time, T / T[0], label = f"$T \\; (\\times{T[0]:.1f} \\; eV)$") -plt.plot(time, Zinj, label = "Injection species") +plt.plot(time, n_e / n_e[0], label = f"$n_e$ ($\\times {n_e[0] * 1e-6:.0e} \\; cm^{{-3}})$") +plt.plot(time, T_i / T_i[0], label = f"$T \\; (\\times{T_i[0]:.1f} \\; eV)$") +plt.plot(time, Zinj, label = "Zinj") +plt.plot(time, Z_Tne, label = "Z_Tne") plt.semilogy() plt.legend() plt.show() diff --git a/scripts_python/plotCumF.py b/scripts_python/plotCumF.py index 3df5eb3..01a8d31 100644 --- a/scripts_python/plotCumF.py +++ b/scripts_python/plotCumF.py @@ -16,7 +16,7 @@ m_i = 1.9712e-25 # paths = ['../2024-12-10_18.45.17/', '../Poisson_50ns_T30Z11/'] # paths = ['../2024-12-11_12.38.27/', '../Poisson_polytropic_fa_T30Z11/', '../Poisson_fa_T30Z11/'] # paths = ['../Poisson_partialAblation/','../Poisson_partialAblation_lowerT/','../Poisson_partialAblation_lowT/','../Poisson_partialAblation_highT/'] -paths = ['../polytropic_80ns_T60/'] +paths = ['../polytropic_80ns_T30/'] labels = [path[3:-1] for path in paths] for path, label in zip(paths, labels): @@ -26,10 +26,12 @@ for path, label in zip(paths, labels): sumF = np.zeros(len(v)) for Z in Zlist: - filesCum_i = sorted(glob.glob(path+'time_*_Z{:.0f}000_fCum_i.csv'.format(Z))) + filename='time_*_Z_{:.1f}_fCum_i.csv'.format(Z) + filesCum_i = sorted(glob.glob(path+filename)) time, x, v, f_i = readF.read(filesCum_i[-1]) sumF += f_i[0] plt.plot(v**2*m_i*0.5/e, f_i[0]*e/m_i/v, label=Z) + print(time) plt.plot(v**2*m_i*0.5/e, sumF*e/m_i/v, label='sum', color='k', linestyle='dashed') plt.yscale('log') diff --git a/scripts_python/plotMom.py b/scripts_python/plotMom.py index 23e3648..2f1dacb 100644 --- a/scripts_python/plotMom.py +++ b/scripts_python/plotMom.py @@ -12,15 +12,15 @@ from scipy.constants import e, k # paths = ['../quasiNeutral_partialAblation/','../Poisson_partialAblation/'] # paths = ['../2024-10-02_14.30.44/'] # paths = ['../quasiNeutral_fullAblation/','../Poisson_fullAblation/'] -paths = ['../2025-04-08_09.36.52/'] +paths = ['../2025-04-09_16.45.52/'] labels = [path[3:-1] for path in paths] for path, label in zip(paths, labels): Zlist = readZlist.read(path+'ZList.csv') filesPhi = sorted(glob.glob(path+'time_*_phi.csv')) - start = 0 - end = 20#len(filesPhi) - every = 5 + start = 80 + end = 85#len(filesPhi) + every = 1 fig, ax = plt.subplots(4, sharex='all') ax[1].set_yscale('log') ax[1].set_ylim(bottom=1e10, top=1e24) @@ -30,7 +30,8 @@ for path, label in zip(paths, labels): ave_ui = np.zeros(len(r)) ave_Ti = np.zeros(len(r)) for Z in Zlist: - filesMom_i = sorted(glob.glob(path+'time_*_Z{:.0f}000_mom_i.csv'.format(Z))) + filename='time_*_Z_{:.1f}_mom_i.csv'.format(Z) + filesMom_i = sorted(glob.glob(path+filename)) fileMom_i = filesMom_i[t] time, r, n_i, u_i, T_i, Zave = readMom.read(fileMom_i) sum_Zni += Zave*n_i diff --git a/scripts_python/readBC.py b/scripts_python/readBC.py index b3de37f..1ed9474 100644 --- a/scripts_python/readBC.py +++ b/scripts_python/readBC.py @@ -4,11 +4,12 @@ def read(filename): # Get time df = pandas.read_csv(filename) time = df['t (s)'].to_numpy() - n = df['n (m^-3)'].to_numpy() - u = df['u (m s^-1)'].to_numpy() - T = df['T (eV)'].to_numpy() + n = df['n_i (m^-3)'].to_numpy() + u = df['u_i (m s^-1)'].to_numpy() + T = df['T_i (eV)'].to_numpy() Zinj = df['Zinj'].to_numpy() + Z_Tne = df['Z_Tne'].to_numpy() - return time, n, u, T, Zinj + return time, n, u, T, Zinj, Z_Tne diff --git a/vlaplex.f90 b/vlaplex.f90 index 1cc9224..4682862 100644 --- a/vlaplex.f90 +++ b/vlaplex.f90 @@ -166,7 +166,7 @@ program VlaPlEx write(*, '(A,ES0.4e3)') 'CFL: ', dt*vf/dr nzMin = 3 - nzMax = 14 + nzMax = 13 nz = nzMax - nzMin + 1 nz = nz + 1 ! Add bin for low Z plasma ! Allocate vectors @@ -189,7 +189,7 @@ program VlaPlEx n_e = 0.0_dp T_e = 0.0_dp Zave = 0.0_dp - Z_list(1) = 1.0_dp ! Low Z bin + Z_list(1) = 0.1_dp ! Low Z bin do iz = nzMin, nzMax Z_list(iz-nzMin+1+1) = float(iz) end do @@ -255,23 +255,17 @@ program VlaPlEx do t = 1, nt time = t * dt + t0 - ! Find new \bar{Z}_i based on T and density + ! Get boundary conditions for specific time call boundaryConditions%get(time, n_bc, u_bc, Temp_bc) - ! Reset previous value - Zave_bc_old = 0.0_dp - ! Initial guess based on average table - call TtoZ%get(Temp_bc, Zave_bc) + ! Find new \bar{Z}_i based on T_e = Temp_bc and n_e = n_bc + call TtoZne%get(Temp_bc, n_bc, Zave_bc) + ! Assign Z(T,n) to bin z_inj = minloc(abs(Z_list - Zave_bc),1) - Zave_bc = Z_list(z_inj) - ! Start iterative process based on T, n_e table - do while (Zave_bc - Zave_bc_old > 0.1_dp) - Zave_bc_old = Zave_bc - call TtoZne%get(Temp_bc, Zave_bc * n_bc, Zave_bc) - end do - z_inj = minloc(abs(Z_list - Zave_bc),1) - Zave_bc = Z_list(z_inj) - u_bc = sqrt(Zave_bc * Temp_bc) - call writeOutputBoundary(t, dt, n_bc, u_bc, Temp_bc, Zave_bc) + ! Calculate inject (sonic) speed + u_bc = sqrt(Z_list(z_inj)* Temp_bc) + ! Calculate ion density to inject + n_bc = n_bc / Z_list(z_inj) + call writeOutputBoundary(t, dt, n_bc, u_bc, Temp_bc, Zave_bc, Z_list(z_inj)) ! f0(j0:nv) = v(j0:nv)**2 / sqrt(PI*Temp_bc**3) * exp(-(v(j0:nv) - u_bc)**2 / Temp_bc) f0(j0:nv) = 1.0_dp / sqrt(PI*Temp_bc) * exp(-(v(j0:nv) - u_bc)**2 / Temp_bc)