From fa59dad8b886c2f15f83a8279a961a491517e60c Mon Sep 17 00:00:00 2001 From: David A Minton Date: Tue, 15 Nov 2022 18:11:40 -0500 Subject: [PATCH] Fixed problems when using EL input --- src/netcdf/netcdf.f90 | 35 ++++++++++++++++++++++++++++------- src/orbel/orbel.f90 | 1 - 2 files changed, 28 insertions(+), 8 deletions(-) diff --git a/src/netcdf/netcdf.f90 b/src/netcdf/netcdf.f90 index 7075d84bc..2cca37386 100644 --- a/src/netcdf/netcdf.f90 +++ b/src/netcdf/netcdf.f90 @@ -619,28 +619,39 @@ module function netcdf_read_frame_system(self, iu, param) result(ierr) end if if ((param%in_form == EL) .or. (param%in_form == XVEL)) then - call check( nf90_get_var(iu%ncid, iu%a_varid, rtemp, start=[1, tslot]), "netcdf_read_frame_system nf90_getvar a_varid" ) + if (.not.allocated(pl%a)) allocate(pl%a(npl)) + if (.not.allocated(tp%a)) allocate(tp%a(ntp)) if (npl > 0) pl%a(:) = pack(rtemp, plmask) if (ntp > 0) tp%a(:) = pack(rtemp, tpmask) call check( nf90_get_var(iu%ncid, iu%e_varid, rtemp, start=[1, tslot]), "netcdf_read_frame_system nf90_getvar e_varid" ) + if (.not.allocated(pl%e)) allocate(pl%e(npl)) + if (.not.allocated(tp%e)) allocate(tp%e(ntp)) if (npl > 0) pl%e(:) = pack(rtemp, plmask) if (ntp > 0) tp%e(:) = pack(rtemp, tpmask) call check( nf90_get_var(iu%ncid, iu%inc_varid, rtemp, start=[1, tslot]), "netcdf_read_frame_system nf90_getvar inc_varid" ) + if (.not.allocated(pl%inc)) allocate(pl%inc(npl)) + if (.not.allocated(tp%inc)) allocate(tp%inc(ntp)) if (npl > 0) pl%inc(:) = pack(rtemp, plmask) if (ntp > 0) tp%inc(:) = pack(rtemp, tpmask) call check( nf90_get_var(iu%ncid, iu%capom_varid, rtemp, start=[1, tslot]), "netcdf_read_frame_system nf90_getvar capom_varid" ) + if (.not.allocated(pl%capom)) allocate(pl%capom(npl)) + if (.not.allocated(tp%capom)) allocate(tp%capom(ntp)) if (npl > 0) pl%capom(:) = pack(rtemp, plmask) if (ntp > 0) tp%capom(:) = pack(rtemp, tpmask) call check( nf90_get_var(iu%ncid, iu%omega_varid, rtemp, start=[1, tslot]), "netcdf_read_frame_system nf90_getvar omega_varid" ) + if (.not.allocated(pl%omega)) allocate(pl%omega(npl)) + if (.not.allocated(tp%omega)) allocate(tp%omega(ntp)) if (npl > 0) pl%omega(:) = pack(rtemp, plmask) if (ntp > 0) tp%omega(:) = pack(rtemp, tpmask) call check( nf90_get_var(iu%ncid, iu%capm_varid, rtemp, start=[1, tslot]), "netcdf_read_frame_system nf90_getvar capm_varid" ) + if (.not.allocated(pl%capm)) allocate(pl%capm(npl)) + if (.not.allocated(tp%capm)) allocate(tp%capm(ntp)) if (npl > 0) pl%capm(:) = pack(rtemp, plmask) if (ntp > 0) tp%capm(:) = pack(rtemp, tpmask) @@ -740,6 +751,10 @@ module function netcdf_read_frame_system(self, iu, param) result(ierr) call self%read_particle_info(iu, param, plmask, tpmask) + if (param%in_form == "EL") then + call pl%el2xv(cb) + call tp%el2xv(cb) + end if ! if this is a GR-enabled run, check to see if we got the pseudovelocities in. Otherwise, we'll need to generate them. if (param%lgr .and. .not.(iu%lpseudo_vel_exists)) then call pl%set_mu(cb) @@ -984,21 +999,24 @@ module subroutine netcdf_read_particle_info_system(self, iu, param, plmask, tpma if (status == nf90_noerr) then call check( nf90_get_var(iu%ncid, iu%origin_xhx_varid, rtemp_arr(1,:)), "netcdf_read_particle_info_system nf90_getvar origin_xhx_varid" ) else - call check( nf90_get_var(iu%ncid, iu%xhx_varid, rtemp_arr(1,:)), "netcdf_read_particle_info_system nf90_getvar xhx_varid" ) + ! [TODO]: This doesn't work when the input mode is EL. This needs to get filled in later after xv2el has been called + ! call check( nf90_get_var(iu%ncid, iu%xhx_varid, rtemp_arr(1,:)), "netcdf_read_particle_info_system nf90_getvar xhx_varid" ) end if status = nf90_inq_varid(iu%ncid, ORIGIN_XHY_VARNAME, iu%origin_xhy_varid) if (status == nf90_noerr) then call check( nf90_get_var(iu%ncid, iu%origin_xhy_varid, rtemp_arr(2,:)), "netcdf_read_particle_info_system nf90_getvar origin_xhy_varid" ) else - call check( nf90_get_var(iu%ncid, iu%xhy_varid, rtemp_arr(2,:)), "netcdf_read_particle_info_system nf90_getvar xhy_varid" ) + ! [TODO]: This doesn't work when the input mode is EL. This needs to get filled in later after xv2el has been called + ! call check( nf90_get_var(iu%ncid, iu%xhy_varid, rtemp_arr(2,:)), "netcdf_read_particle_info_system nf90_getvar xhy_varid" ) end if status = nf90_inq_varid(iu%ncid, ORIGIN_XHZ_VARNAME, iu%origin_xhz_varid) if (status == nf90_noerr) then call check( nf90_get_var(iu%ncid, iu%origin_xhz_varid, rtemp_arr(3,:)), "netcdf_read_particle_info_system nf90_getvar origin_xhz_varid" ) else - call check( nf90_get_var(iu%ncid, iu%xhz_varid, rtemp_arr(3,:)), "netcdf_read_particle_info_system nf90_getvar xhz_varid" ) + ! [TODO]: This doesn't work when the input mode is EL. This needs to get filled in later after xv2el has been called + ! call check( nf90_get_var(iu%ncid, iu%xhz_varid, rtemp_arr(3,:)), "netcdf_read_particle_info_system nf90_getvar xhz_varid" ) end if do i = 1, npl @@ -1012,21 +1030,24 @@ module subroutine netcdf_read_particle_info_system(self, iu, param, plmask, tpma if (status == nf90_noerr) then call check( nf90_get_var(iu%ncid, iu%origin_vhx_varid, rtemp_arr(1,:)), "netcdf_read_particle_info_system nf90_getvar origin_vhx_varid" ) else - call check( nf90_get_var(iu%ncid, iu%vhx_varid, rtemp_arr(1,:)), "netcdf_read_particle_info_system nf90_getvar vhx_varid" ) + ! [TODO]: This doesn't work when the input mode is EL. This needs to get filled in later after xv2el has been called + ! call check( nf90_get_var(iu%ncid, iu%vhx_varid, rtemp_arr(1,:)), "netcdf_read_particle_info_system nf90_getvar vhx_varid" ) end if status = nf90_inq_varid(iu%ncid, ORIGIN_VHY_VARNAME, iu%origin_vhy_varid) if (status == nf90_noerr) then call check( nf90_get_var(iu%ncid, iu%origin_vhy_varid, rtemp_arr(2,:)), "netcdf_read_particle_info_system nf90_getvar origin_vhy_varid" ) else - call check( nf90_get_var(iu%ncid, iu%vhy_varid, rtemp_arr(2,:)), "netcdf_read_particle_info_system nf90_getvar vhy_varid" ) + ! [TODO]: This doesn't work when the input mode is EL. This needs to get filled in later after xv2el has been called + ! call check( nf90_get_var(iu%ncid, iu%vhy_varid, rtemp_arr(2,:)), "netcdf_read_particle_info_system nf90_getvar vhy_varid" ) end if status = nf90_inq_varid(iu%ncid, ORIGIN_VHZ_VARNAME, iu%origin_vhz_varid) if (status == nf90_noerr) then call check( nf90_get_var(iu%ncid, iu%origin_vhz_varid, rtemp_arr(3,:)), "netcdf_read_particle_info_system nf90_getvar origin_vhz_varid" ) else - call check( nf90_get_var(iu%ncid, iu%vhz_varid, rtemp_arr(3,:)), "netcdf_read_particle_info_system nf90_getvar vhz_varid" ) + ! [TODO]: This doesn't work when the input mode is EL. This needs to get filled in later after xv2el has been called + ! call check( nf90_get_var(iu%ncid, iu%vhz_varid, rtemp_arr(3,:)), "netcdf_read_particle_info_system nf90_getvar vhz_varid" ) end if do i = 1, npl diff --git a/src/orbel/orbel.f90 b/src/orbel/orbel.f90 index 9dc68c397..5e7c4a989 100644 --- a/src/orbel/orbel.f90 +++ b/src/orbel/orbel.f90 @@ -30,7 +30,6 @@ module subroutine orbel_el2xv_vec(self, cb) call orbel_el2xv(self%mu(i), self%a(i), self%e(i), self%inc(i), self%capom(i), & self%omega(i), self%capm(i), self%xh(:, i), self%vh(:, i)) end do - deallocate(self%a, self%e, self%inc, self%capom, self%omega, self%capm) return end subroutine orbel_el2xv_vec