research/incomp-flame-1d
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

added reaction rate function

Browse source
This commit is contained in:
ignis 2019-08-30 07:49:27 +09:00
parent 042e594e07
commit b7fc2451a3
1 changed files with 37 additions and 7 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

44
code/ysolve.f90
View file

@ -380,13 +380,7 @@
CALL d2fnonp(nx,hx,ux(:,:),d2ux(:,:),2,1)
DO i=1,nx
wrate=pre*ux(1,i)*exp(-ac/(1.+bc*(ux(2,i)))) !wrate
IF ((ux(2,i)).le.c_ref) THEN
wrate=min_wr
IF ((ux(2,i)).gt.c_cut) wrate=((refwr-min_wr)*exp(prof_wr*(ux(2,i)-c_ref))+ &
min_wr-refwr*exp(prof_wr*(c_cut-c_ref)))/(1.-exp(prof_wr*(c_cut-c_ref)))
ENDIF
wrate=rate_1step(ux(1,i), ux(2,i))
! -0.5*( d(u*Y)/dx + u*dY/dx + Y*du/dx ) + D*d2Y/d2x
f(i,1) = - ( u(i)*dux(1,i) ) + diff * d2ux(1,i) - wrate
@ -416,5 +410,41 @@
f(nx,2)=f(nx,2)-Dt
END SUBROUTINE fns
!------------------------------------------------------------------------
real function rate_1step (yr, theta)
real, intent(in) :: yr
real, intent(in) :: theta
real :: y
real :: t_reduce
y=yr
if(yr.lt.0.) y=0.
if(yr.gt.1.) y=1.
t_reduce=theta
if(theta.lt.0.) t_reduce=0.
if(theta.gt.1.) t_reduce=1.
if (t_reduce.gt.c_ref) then
rate_1step = pre*y*exp(-ac/(1.+bc*t_reduce))
else if (t_reduce.le.c_cut) then
rate_1step = min_wr
else
rate_1step = &
((refwr-min_wr)*exp(prof_wr*(t_reduce-c_ref)) + min_wr - refwr*exp(prof_wr*(c_cut-c_ref))) &
/ (1.-exp(prof_wr*(c_cut-c_ref)))
endif
end function rate_1step
END MODULE ysolve