move initial condition selecting to upper level

This commit is contained in:
ignis 2020-02-26 13:15:57 +09:00
parent 4c4411d6d3
commit b68cadc136

View file

@ -90,22 +90,25 @@
SUBROUTINE solve
IF (reaction_type == "onestep") THEN
CALL solve_ (fonestep)
CALL solve_ (fonestep, SET_IC_ONESTEP)
ELSE IF (reaction_type == "twostep") THEN
CALL solve_ (fns)
CALL solve_ (fns, SET_IC_TWOSTEP)
ELSE
WRITE(*,*)'ERROR, UNSUPPORTED REACTION ', trim(reaction_type), '. STOPPING'
STOP
END IF
END SUBROUTINE solve
SUBROUTINE solve_ (rhs)
SUBROUTINE solve_ (rhs, SET_IC)
INTERFACE
SUBROUTINE rhs(r1,f)
REAL, INTENT(IN), DIMENSION(:,:) :: r1
REAL, INTENT(OUT), DIMENSION(:,:) :: f
END SUBROUTINE rhs
SUBROUTINE SET_IC
END SUBROUTINE SET_IC
END INTERFACE
INTEGER :: i,j,k,savenum
@ -114,6 +117,12 @@
CALL SET_IC
pflame=0.
DO i=1,nx
pflame=pflame+y1(i,fctrl_species)*hx
ENDDO
t_now=0.
t_uf=0.
@ -331,7 +340,7 @@
END DO
END SUBROUTINE SET_BC
SUBROUTINE SET_IC
SUBROUTINE BASE_FLAME_PROFILE(x)
INTEGER :: i, ifl, si
REAL :: xi
REAL :: x(nx)
@ -355,18 +364,39 @@
ELSE
x(i)=0.5+REAL(i-nx+si)/REAL(ifl)
ENDIF
END DO
if ( reaction_type == "onestep" ) then
! initialize species field
! max_ysum = maxval(y1(:,1) + y1(:,2))
! y1(:,1) = y1(:,1) / maxval(y1(:,1) + y1(:,2))
! y1(:,2) = y1(:,2) / maxval(y1(:,1) + y1(;,2))
! y1(:,1) = y1(:,1) / (y1(1,1) + y1(1,2))
! y1(:,2) = y1(:,2) / (y1(1,1) + y1(1,2))
END SUBROUTINE BASE_FLAME_PROFILE
SUBROUTINE SET_IC_ONESTEP
INTEGER :: i
REAL :: xi
REAL :: x(nx)
CALL BASE_FLAME_PROFILE(x)
DO i=1,nx
xi = x(i)
y1(i,1)=(1.-xi) ! reactant mass fraction
y1(i,2)=(xi) ! reactant mass fraction
ENDDO
else if ( reaction_type == "twostep" ) then
END SUBROUTINE SET_IC_ONESTEP
SUBROUTINE SET_IC_TWOSTEP
INTEGER :: i
REAL :: xi
REAL :: x(nx)
CALL BASE_FLAME_PROFILE(x)
DO i=1,nx
xi = x(i)
@ -377,26 +407,7 @@
y1(i,3) = xi
ENDDO
else
WRITE(*,*) 'ERROR, UNDEFINED REACTION TYPE ', reaction_type
stop
end if
! initialize species field
! max_ysum = maxval(y1(:,1) + y1(:,2))
! y1(:,1) = y1(:,1) / maxval(y1(:,1) + y1(:,2))
! y1(:,2) = y1(:,2) / maxval(y1(:,1) + y1(;,2))
! y1(:,1) = y1(:,1) / (y1(1,1) + y1(1,2))
! y1(:,2) = y1(:,2) / (y1(1,1) + y1(1,2))
! calculate flame position
pflame=0.
DO i=1,nx
pflame=pflame+y1(i,fctrl_species)*hx
ENDDO
END SUBROUTINE SET_IC
END SUBROUTINE SET_IC_TWOSTEP
SUBROUTINE init_solver