Time derivative of Euler field, the residuals function.
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(euler_1D_openmp), | intent(in) | :: | self | Euler field. |
||
real(kind=R_P), | intent(in), | optional | :: | t | Time. |
Euler field time derivative.
function dEuler_dt(self, t) result(dState_dt)
!---------------------------------------------------------------------------------------------------------------------------------
!< Time derivative of Euler field, the residuals function.
!---------------------------------------------------------------------------------------------------------------------------------
class(euler_1D_openmp), intent(IN) :: self !< Euler field.
real(R_P), optional, intent(IN) :: t !< Time.
class(integrand), allocatable :: dState_dt !< Euler field time derivative.
real(R_P), allocatable :: F(:,:) !< Fluxes of conservative variables.
real(R_P), allocatable :: P(:,:) !< Primitive variables.
real(R_P), allocatable :: PR(:,:,:) !< Left (1) and right (2) (reconstructed) interface values of primitive variables.
integer(I_P) :: i !< Counter.
!---------------------------------------------------------------------------------------------------------------------------------
!---------------------------------------------------------------------------------------------------------------------------------
allocate(F(1:self%Nc, 0:self%Ni))
!$OMP PARALLEL DO PRIVATE(i) SHARED(self, F)
do i=0, self%Ni
F(:, i) = 0._R_P
enddo
allocate(P(1:self%Np, 1-self%Ng:self%Ni+self%Ng))
!$OMP PARALLEL DO PRIVATE(i) SHARED(self, P)
do i=1-self%Ng, self%Ni+self%Ng
P(:, i) = 0._R_P
enddo
allocate(PR(1:self%Np, 1:2, 0:self%Ni+1))
!$OMP PARALLEL DO PRIVATE(i) SHARED(self, P)
do i=0, self%Ni+1
PR(:, :, i) = 0._R_P
enddo
! compute primitive variables
!$OMP PARALLEL DO PRIVATE(i) SHARED(self, P)
do i=1, self%Ni
P(:, i) = self%conservative2primitive(self%U(:, i))
enddo
call self%impose_boundary_conditions(primitive=P)
call self%reconstruct_interfaces_states(primitive=P, r_primitive=PR)
! compute fluxes by solving Rimeann Problems at each interface
!$OMP PARALLEL DO PRIVATE(i) SHARED(self, F, PR)
do i=0, self%Ni
call self%riemann_solver(r1=PR(self%Ns+3, 2, i ), &
u1=PR(self%Ns+1, 2, i ), &
p1=PR(self%Ns+2, 2, i ), &
g1=PR(self%Ns+4, 2, i ), &
r4=PR(self%Ns+3, 1, i+1), &
u4=PR(self%Ns+1, 1, i+1), &
p4=PR(self%Ns+2, 1, i+1), &
g4=PR(self%Ns+4, 1, i+1), &
F=F(:, i))
if (self%Ns>1) then
if (F(1, i)>0._R_P) then
F(1:self%Ns, i) = PR(1:self%Ns, 2, i )/PR(self%Ns+3, 2, i )*F(1, i)
else
F(1:self%Ns, i) = PR(1:self%Ns, 1, i+1)/PR(self%Ns+3, 1, i+1)*F(1, i)
endif
endif
enddo
! compute residuals
allocate(euler_1D_openmp :: dState_dt)
select type(dState_dt)
class is(euler_1D_openmp)
dState_dt = self
endselect
!$OMP PARALLEL PRIVATE(i) SHARED(self, dState_dt, F)
select type(dState_dt)
class is(euler_1D_openmp)
!$OMP DO
do i=1, self%Ni
dState_dt%U(:, i) = (F(:, i - 1) - F(:, i)) / self%Dx
enddo
endselect
!$OMP END PARALLEL
return
!---------------------------------------------------------------------------------------------------------------------------------
endfunction dEuler_dt