Boundary Conditions and Initial Conditions

This commit is contained in:
ignis 2020-02-21 08:04:25 +09:00
parent 9e66226c65
commit 964f8ccda6
2 changed files with 65 additions and 17 deletions

View file

@ -14,6 +14,18 @@ real, private :: hrp_twostep
contains
subroutine init_chemistry
character(len=40) :: nrxn_string
if (nrxn == 1) then
reaction_type = "onestep"
else if (nrxn == 2) then
reaction_type = "twostep"
else
write(nrxn_string, *) nrxn
reaction_type = trim(nrxn_string) // "-step"
end if
if ( reaction_type == "onestep" ) then
lambda_onestep = pre * exp ( - beta / hrp )
else if ( reaction_type == "twostep" ) then
@ -22,6 +34,7 @@ subroutine init_chemistry
beta1_twostep = beta1
hrp_twostep = hrp
else
WRITE(*,*) 'ERROR, UNDEFINED REACTION TYPE ', reaction_type
stop
end if

View file

@ -6,6 +6,7 @@
IMPLICIT NONE
REAL, DIMENSION(:), ALLOCATABLE :: u
REAL, DIMENSION(:), ALLOCATABLE :: inletbc
REAL, DIMENSION(:,:), ALLOCATABLE :: y1,y2,yf,yold
REAL, DIMENSION(:), ALLOCATABLE :: uxt,duxt
REAL, DIMENSION(:), ALLOCATABLE :: dm
@ -121,7 +122,7 @@
pflame=0.
DO i=1,nx
pflame=pflame+y1(i,1)*hx
uxt(i)=y1(i,3)
uxt(i)= 1. - y1(i,1)
ENDDO
CALL dfnonp(nx,hx,uxt,duxt,1,1)
@ -287,14 +288,17 @@
END SUBROUTINE write_pre
SUBROUTINE SET_BC
y1(1,1)=1.
y1(1,2)=0.
y1(1,3)=0.
INTEGER :: i
DO i = 1, nsp
y1(1,i) = inletbc(i)
END DO
END SUBROUTINE SET_BC
SUBROUTINE SET_IC
INTEGER :: i, ifl, si
REAL :: xi
REAL :: x(nx)
REAL :: max_ysum
! initial flame thickness (0.2 pi n grids)
@ -309,22 +313,37 @@
! initialize Yr field
DO i=1,nx
IF(i< nx-(si+ifl/2)) THEN
xi=0.+ctmp
x(i)=0.+ctmp
ELSE IF(i> nx-(si-ifl/2)) THEN
xi=1.
x(i)=1.
ELSE
xi=0.5+REAL(i-nx+si)/REAL(ifl)
x(i)=0.5+REAL(i-nx+si)/REAL(ifl)
ENDIF
! y1(i,1,1)=(1.-xi)*1.
y1(i,1)=(1.-xi) ! reactant mass fraction
END DO
y1(i,2) = (1./2.) * (lambda1/lambda2) * y1(i,1) * &
exp (-(beta1*(1. - xi))/(1. - hrp*(1. - xi)))
if ( reaction_type == "onestep" ) then
DO i=1,nx
xi = x(i)
y1(i,3) = xi
y1(i,1)=(1.-xi) ! reactant mass fraction
ENDDO
else if ( reaction_type == "twostep" ) then
DO i=1,nx
xi = x(i)
y1(i,1)=(1.-xi) ! reactant mass fraction
y1(i,2) = (1./2.) * (lambda1/lambda2) * y1(i,1) * &
exp (-(beta1*(1. - xi))/(1. - hrp*(1. - xi)))
y1(i,3) = xi
ENDDO
else
WRITE(*,*) 'ERROR, UNDEFINED REACTION TYPE ', reaction_type
stop
end if
ENDDO
! initialize species field
! max_ysum = maxval(y1(:,1) + y1(:,2))
@ -350,14 +369,29 @@
CALL init_chemistry
ALLOCATE( u(nx),STAT=ierr) ; u=0.
ALLOCATE( y1(nx,3),STAT=ierr) ; y1=0.
ALLOCATE( y2(nx,3),STAT=ierr) ; y2=0.
ALLOCATE( yf(nx,3),STAT=ierr) ; yf=0.
ALLOCATE(yold(nx,3),STAT=ierr) ; yold=0.
ALLOCATE( inletbc(nsp),STAT=ierr) ; inletbc=0.
ALLOCATE( y1(nx,nsp),STAT=ierr) ; y1=0.
ALLOCATE( y2(nx,nsp),STAT=ierr) ; y2=0.
ALLOCATE( yf(nx,nsp),STAT=ierr) ; yf=0.
ALLOCATE(yold(nx,nsp),STAT=ierr) ; yold=0.
ALLOCATE(uxt(nx),STAT=ierr) ; uxt=0.
ALLOCATE(duxt(nx),STAT=ierr) ; duxt=0.
ALLOCATE(dm(nx),STAT=ierr) ; dm=diff
if ( reaction_type == "onestep" ) then
inletbc(1) = 1.
else if ( reaction_type == "twostep" ) then
inletbc(1) = 1.
inletbc(2) = 0.
inletbc(3) = 0.
else
WRITE(*,*) 'ERROR, UNDEFINED REACTION TYPE ', reaction_type
stop
end if
CALL ludcmp(nx,100,100,1,0,0)
END SUBROUTINE init_solver
@ -365,6 +399,7 @@
SUBROUTINE finalize_solver
DEALLOCATE( u)
DEALLOCATE(inletbc)
DEALLOCATE(y1)
DEALLOCATE(y2)
DEALLOCATE(yf)