program EDDA
	! EDDA:integrated simulation of debris flow erosion, deposition and property changes
	! by H.X. Chen, HKUST, 2015
	use input_file_defs; use input_vars
	use grids; use model_vars;use output_file_defs
	implicit none
	integer::grd
	integer:: i,ii,j,k,l,m,n,imx1,mnd,maxzo,maxdirection
	integer:: nodata,sctr
	integer:: ncol,nrow,u(25)
	character (len=255):: outfil 
	character (len=14):: header(6)
	character (len=8):: date
	character (len=10):: time
	character (len=4):: stp
	character (len=31):: scratch,irfil
	character (len=7):: vrsn
	character (len=11):: bldate


subroutine dfs(imx1,nrow,ncol,header,u,maxdirection)
	! subroutine for debris flow simulation
	use model_vars;use input_vars; use grids;use input_file_defs;use output_file_defs
	use inflow_vars
	implicit none
	integer:: i,ii,j,k,m,n,nsteps,ncol,nrow,trace,imx1,nt,nq,inflowsize,vv,vvmax,maxdirection
	integer:: u(25)
	double precision:: tnext
	double precision:: tempri(imx1),tempinflowh(imx1),tempinflowrho(imx1)
	double precision:: fvpredi(imx1,maxdirection),fvpredi2(imx1,maxdirection),fhpredi(imx1),frhopredi(imx1),cv(imx1),erorate(imx1),deporate(imx1),fhpredi2(imx1),frhopredi2(imx1),fybar(imx1,maxdirection),fq(imx1,maxdirection)
	double precision:: sf,sfy,sfmiu,sfmanning,coemiu,coemanning,tao,taoc,grav,normfric,absubar(imx1),normfriccoe,miudebris,gammadeb,rhoero,rhodepo(imx1),manningb,manningm
	double precision:: tol
	double precision:: hi,hn,grad,ybar,agrad,cvbar,normfriccoebar,miubar,hbar,manningbar,frhobar,frhoflux
	double precision:: parai,paran,aa,bb,cc,dd,vel,qq(imx1,maxdirection),qnet,hinflow,cellarea,qqmass(imx1,maxdirection),qmassnet,dttest,dfhtest,dpfhtest,fvlimit,dt0,localvdiff,dv
	double precision:: hour,cosslope
	character (len=25):: stp
	character (len=255):: outfil
	character (len=14)::  header(6)
	double precision:: totalvolume,totalflowvolume,totaldepositvolume
	double precision:: lambdainverse(imx1),fvdepo(imx1),tanthetae,sinthetae,cvlimit(imx1)
	double precision:: tempinierodithick(imx1),tempdebdepothick(imx1)
	double precision:: artivis
	double precision:: width,rholimit(imx1),vx,vy,fhtemp,fvtemp,cvtemp


subroutine flodir(grd,pth,ncol,nrow,celsiz,nodat,pf,pf1,ctr,imax,temp,u,infil,param,header,u1,fp,outflow)
	implicit none
	integer grd,pth,i,m,ncol,nrow,ctr,imax,u,u1
	double precision param(6),nodat,celsiz,cns,cew,temp(pth)
	double precision east,west,north,south
	double precision pf(imax),pf1(grd),nodats
	integer fp(imax,8),cellrecord(grd),cellnow
	character*14 header(6)
	character*255 infil
	logical outflow(imax)
	
	
module grids
	integer,allocatable:: pf2(:),sctrtrace(:),indx(:),nxt(:),nv(:),nvu(:),dir(:),road(:),catchment(:),mouthpoint(:),slopedeposittrace(:),debrisflowtrace(:),warninglevel(:),warninglevelave(:)
	integer,allocatable:: dsctr(:),dsc(:),zo(:),temp(:),group(:),maxelecell(:),finalcell(:)
	double precision,allocatable:: rikzero(:),stemp(:)
	double precision,allocatable::  rik(:),rik1(:),ri(:),rizero(:),pf1(:),pfmin(:)
	double precision,allocatable:: ro(:),wf(:),ir(:),tfg(:)
	double precision,allocatable:: zmax(:),slo(:),depth(:),ele(:),finalvolume(:),finaldistance(:),finaldeltah(:),ltstar(:),debrisdeposit(:)
	double precision,allocatable:: zfmin(:),fsmin(:),pmin(:),fdepth(:),maxele(:),volume(:),slopedeposit(:),sloperoadrisk(:),cumurain(:),debrisroadrisk(:)
	double precision,allocatable:: maxindirisk(:),aveindirisk(:)
	double precision,allocatable:: pfmax(:,:)
	
	! parameters for whole processes of debris flow simulation
	integer, allocatable:: fp(:,:)
	double precision, allocatable:: fv(:,:),maxfv(:)
	double precision, allocatable:: rideb(:,:)
	double precision, allocatable:: inierodithick(:),debdepothick(:),fh(:),frho(:),maxfh(:)
	logical,allocatable:: outflow(:)
	

subroutine inflow_read(ulog,uinflow)
	use input_file_defs;use inflow_vars
	implicit none
	integer i,j,linct,nopulse
	integer,intent(in)::ulog,uinflow
	logical ans


module inflow_vars
	! parameters for whole processes of debris flow simulation
	integer noinflow
	integer,allocatable:: inflowcell(:)
	double precision inflowperiod,inflowdt
	double precision,allocatable:: inflowht(:,:),inflowhq(:,:),inflowhcv(:,:)


module input_file_defs
  	character (len=255):: heading,slofil,zonfil,zfil,demfil,ltstarfil,roadfil,catchmentfil,mouthpointfil
  	character (len=255):: depfil
  	character (len=255):: rizerofil,nxtfil,ndxfil,dirfil
	character (len=255):: wffil,dscfil,init,title
  	character (len=255), allocatable:: rifil(:)
 	character (len=224):: folder
	character (len=8):: suffix
	
	! files for debris flow simulation
	character (len=255)::manningfil
	character (len=255)::inflowfil


module input_vars
	logical:: ans,outp(8),rodoc,lskip,lany,llus,lps0,unsat0,bkgrof,depositsave,pfsave,risksave,warninglevelsave,tracesave
	logical,allocatable:: unsat(:) 
 	integer:: imax,row,col,nwf,tx,nmax
	integer:: flag,nper,aif
  	integer:: nzsb,nzst,mmax
	integer:: nzon,nout,nsection
	integer,allocatable:: ksav(:)
	double precision :: uww,zmin,t,dep,czmax,crizero,slomin
	double precision,allocatable:: cb(:),sdcb(:),ct(:),sdct(:),phib(:),sdphib(:),phit(:),sdphit(:),phibb(:),sdphibb(:),phibt(:),sdphibt(:),uwsb(:),uwst(:),difb(:),dift(:),ksb(:),kst(:),thsatb(:),thresib(:),thsatt(:),thresit(:),alphat(:),alphab(:),pt(:),pb(:)
	double precision :: rlb,rub,rtb,fl,fu,dx, cltstar, lbstar,qa
	double precision,allocatable:: capt(:),cri(:),tsav(:)
	character (len=5):: flowdir
	
	!parameters for debris flow simulation
	double precision:: alpha1,beta1,alpha2,beta2,kresis,cmanning,limitfr
	double precision,allocatable::manning(:),kero(:)
	double precision:: d50,cvstar,coedepo,cs
	double precision:: dtmin, dtmax,dti,dtd,simul,tout,toldh,toldhp
	logical:: rainsimul,erosionsimul,inflowsimul,sepdepositionsimul
    logical:: flowdepthsave,maxflowdepthsave,fvsave,maxfvsave,fulldyna,erodepthsave,debdepodepthsave,cvsave,totaldepthsave
	
	
subroutine irdgrd(grd,pth,ncol,nrow,celsiz,inodat,mnd,y,y1,y2,ctr,imax,itemp,u,infil,param,header,u1)
	implicit none
	integer grd,pth,i,m,ncol,nrow,ctr,imax,u,inodat,u1,count
	integer y(imax),y1(grd),y2(grd),itemp(pth),mnd
	double precision param(6),celsiz,cns,cew,nodat
	double precision east,west,north,south
	character*14 header(6)
	character*255 infil
	
	
module model_vars
    ! variables for rainfall-induced slope failures
	integer:: nts,kper,nmax1,nmax2,nmn,nmin,maxcatchment,cellcount(100)
	integer,allocatable:: jsav(:),templevel(:),roadindx(:),roadsectioncell(:)
	double precision:: test1,dg2rad,qb,level1,level2,level3,level4,level5
	double precision,allocatable:: q(:)
	double precision:: eps,tmin,tmax,ts,tinc,qt,tns,beta,qmax !
	double precision:: test,nodat,sumex,dusz,dcf,vf0 !(:)
	double precision:: celsiz,param(6),parami(6),ti,tis,pi,smt,lard
	double precision,allocatable:: p(:),ptran(:),pzero(:),bline(:)
	double precision,allocatable:: r(:),fc(:),fw(:),inithzt(:),thzt(:),thzb(:),kz(:),tcap(:),inidesat(:),desat(:)
	double precision,allocatable:: trz(:),uwsp(:),uwspt(:),uwspb(:),gst(:),gsb(:),gs(:),qtime(:)
	double precision,allocatable:: catchmentarea(:),catchmentdeltah(:),catchmentmouth(:)
	
	! variables for debris flow simulation
	
	! other variables
	double precision:: rhow,rhos
	
	! Auxiliary parameters
	double precision:: t0,dt,ttout,tnow,tminimum,tmaximum,tmean
	double precision, allocatable:: dmax(:),tdmax(:)
	integer:: ntsdeb,maxnts,avedt


module output_file_defs
	
	! files for debris flow simulation
	character (len=25)::fhfil='Flow_depth_'
	character (len=25)::maxfvfil='Max_flow_velocity_'
    character (len=25)::fvfil='Flow_velocity_'
	character (len=25)::maxfhfil='Max_flow_depth_'
	character (len=25)::erodepthfil='Erosion_depth_'
	character (len=25)::debdepodepthfil='Deposit_depth_'
	character (len=25)::totaldepthfil='Total_depth_'
	character (len=25)::cvfil='Volumetric_sediment_concentration_'


subroutine srdgrd(grd,pth,ncol,nrow,celsiz,nodat,pf,pf1,ctr,imax,temp,u,infil,param,header,u1)
	implicit none
	integer grd,pth,i,m,ncol,nrow,ctr,imax,u,u1
	double precision param(6),nodat,celsiz,cns,cew,temp(pth)
	double precision east,west,north,south
	double precision pf(imax),pf1(grd),nodats
	character*14 header(6)
	character*255 infil


subroutine ssvgrd(z3,grd,z2,nrow,ncol,u,nodata,param,u1,outfil,ti,header)
	implicit none
	integer i,j,nrow,ncol,u,grd,ctr,u1,ia,m
 !  Map 2-d array, z2, onto 1-d array, z3.
 !  Reverse rows and columns, because Fortran stores arrays
 !  by column.  z2 includes nodata values and is used to map
 !  z3 to the grid format
	double precision z2(ncol,nrow),z3(grd),nodata
	double precision test, param(6),ti,a
	character*14 header(6)
	character*1 sp
	character*255 outfil
 	character*31 scratch
 
 
subroutine trini(ulog,uini,dg2rad)
	! reads initialization file, R.L. Baum, USGS, revised by H.X. Chen, 2014
	use input_file_defs; use input_vars
	implicit none
	integer:: i,j,iz,linct
	integer, intent(in):: ulog,uini
	double precision, intent(in):: dg2rad
	character (len=31):: scratch