JGonzalez/vlaplex
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

Add Z injection index

Browse source
This commit is contained in:
JHendrikx 2025-02-04 17:23:19 +01:00
commit 484e6bea72
1 changed files with 26 additions and 11 deletions
Download patch file Download diff file Expand all files Collapse all files

37
vlaplex.f90
View file

@ -47,7 +47,7 @@ program VlaPlEx
real(dp):: time real(dp):: time
real(dp):: dr, dv, dt real(dp):: dr, dv, dt
integer:: nr, nv, nt, nz integer:: nr, nv, nt, nz
integer:: i, iz, j, t integer:: i, iz, j, t, z_inj
integer:: j0 ! First integer of positive velocity integer:: j0 ! First integer of positive velocity
real(dp):: Temp_bc ! Temperature real(dp):: Temp_bc ! Temperature
@ -241,18 +241,33 @@ program VlaPlEx
time = t * dt + t0 time = t * dt + t0
call boundaryConditions%get(time, n_bc, u_bc, Temp_bc, Zave_bc) call boundaryConditions%get(time, n_bc, u_bc, Temp_bc, Zave_bc)
call writeOutputBoundary(t, dt, n_bc, u_bc, Temp_bc, Zave_bc) call writeOutputBoundary(t, dt, n_bc, u_bc, Temp_bc, Zave_bc)
z_inj = minloc(abs(Z_list - T_to_Z(Temp_bc)),1)
Zave_bc = Z_list(z_inj)
u_bc = sqrt(Zave_bc * Temp_bc) u_bc = sqrt(Zave_bc * Temp_bc)
do iz = 1, nz
! 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)
f0 = f0 * n_bc / (sum(f0)*dv)
! Boundary conditions ! f0(j0:nv) = v(j0:nv)**2 / sqrt(PI*Temp_bc**3) * exp(-(v(j0:nv) - u_bc)**2 / Temp_bc)
! r = r0, v>0 f0(j0:nv) = 1.0_dp / sqrt(PI*Temp_bc) * exp(-(v(j0:nv) - u_bc)**2 / Temp_bc)
f_i_old(iz,1,j0:nv) = f0 f0 = f0 * n_bc / (sum(f0)*dv)
f_i(iz,1,j0:nv) = f_i_old(iz,1,j0:nv)
T_i(iz,1) = Temp_bc !f_i_old(:,1,j0:nv) = 0
end do !f_i(:,1,j0:nv) = 0
!T_i(:,1) = 0
f_i_old(z_inj,1,j0:nv) = f0
f_i(z_inj,1,j0:nv) = f_i_old(z_inj,1,j0:nv)
T_i(z_inj,1) = Temp_bc
!do iz = 1, nz
! ! 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)
! f0 = f0 * n_bc / (sum(f0)*dv)
! ! Boundary conditions
! ! r = r0, v>0
! f_i_old(iz,1,j0:nv) = f0
! f_i(iz,1,j0:nv) = f_i_old(iz,1,j0:nv)
! T_i(iz,1) = Temp_bc
!end do
T_e = Temp_bc T_e = Temp_bc
print *, 'Time: ', time * t_ref print *, 'Time: ', time * t_ref
print *, 'Temp_bc: ', Temp_bc print *, 'Temp_bc: ', Temp_bc