## Boundary condition implementation in Hydrus-1D

A discussion forum for Hydrus-1D users.
Hue
Posts: 6
Joined: Tue Jul 16, 2019 12:32 pm
Location: Germany

### Boundary condition implementation in Hydrus-1D

Dear Hydrus team and Forum Members,

We simulated a simple conservative tracer experiment in a homogeneous soil column (L=100cm), with constant hydraulic head boundary conditions inducing a saturated column with flow from top to bottom.

We are simulating a constant injection of a conservative tracer on the top. We are interested in the break-through curves along the soil profile. We would like to compare the Hydrus output to the Ogata-Banks solution, expecting they should be identical for the entire column length.

We played with the solute boundary conditions and have two cases:
- Case 1, where we apply a constant concentration boundary condition in Hydrus-1D
- Case 2, where we apply a constant flux boundary condition

The breakthrough curves are indeed identical but only at the distance from 0cm to approximately 92cm from the source of injection (case 1). Starting from 92cm till the column outlet, or the final 8cm or so of the profile, the Hydrus solutions and OgataBanks solutions are increasingly deviating (case 1). We have tried several different sytem settings and conditions, and still got the same kind of differences between the numerical and analytical solutions.

Case 1
Hydrus-OgataBanks_1.png (27.59 KiB) Viewed 2300 times

But when the upper BC for the solute transport is concentration flux BC (case 2), the solutions are deviating along the entire column length but gradually converging to the Ogata-Banks solutions at the column outlet. So the Hydrus breakthrough curve at the outlet looks more similar to the analytical solutions than in the constant concentration BC case. We all know that it’s the constant concentration BC in Ogata-Banks.

Case 2
Hydrus-OgataBanks_2.png (41.06 KiB) Viewed 2299 times

We suspect what might have caused these differences in Hydrus solution is in the way how Hydrus implements the boundary conditions which we are not sure of. We would appreciate any comments of clarification and suggestions for minimizing the difference between the analytical- and the Hydrus- solution.

A summary of our model is below, and I also attached the Hydrus codes in the zip file.
L = 100 cm
i = 0.01
K_sat = 8.64 cm/day
Porosity = 0.25
Dispersivity = 4.423 cm
Diffusion coefficient = 0
Simulation time = 600 days
Water flow upper and lower BCs: Constant pressure head
Solute transport upper BC: Constant concentration of 1 mmol/cm3 (case 1); Concentration flux of 1 mmol/cm3 (case 2)
Solute transport lower BC: Zero concentration gradient

Thank you.
Hue
Attachments
cases1-2.zip

Jirka
Posts: 5268
Joined: Sat Mar 16, 2002 3:47 pm
Location: USA
Location: Riverside, CA

### Re: Boundary condition implementation in Hydrus-1D

The solute transport part of HYDRUS-1D has been extensively and successfully tested against the analytical solutions from CFitm, CFitIm, and CXTFIT. Note that these analytical solutions are valid for different surface (Dirichlet or Cauchy) BCs, as well as for different bottom boundary conditions (finite and infinite columns). Obviously, the solution must to be different for Dirichlet and Cauchy BCs (as more solute enters the column for the Dirichlet BC), and for different bottom BCs. You show only one analytical solution, which you call the Ogata-Banks solution (which I’m not familiar with), without actually specifying for which BCs it was derived. You clearly cannot compare two different numerical solutions for different surface BCs with a single analytical solution. J.

Hue
Posts: 6
Joined: Tue Jul 16, 2019 12:32 pm
Location: Germany

### Re: Boundary condition implementation in Hydrus-1D

Thank you, Jirka, for the quick comment. We will be looking at the 3 analytical solutions you mentioned.

With our cases, I should have been clearer.
As an analytical solution to the advection-dispersion equation (ADE), the Ogata-Banks solution (1961) we used has the following form:
C(x,t) = (Co/2) * [erfc[(x-vt)/2sqrt(Dt)] + exp(vx/D) * erfc[(x+vt)/2sqrt(Dt)]]
Where erfc is the complementary error function, and:
C = solute concentration [M/L^3]
x = distance of interest from the source of solute injection, along the flow path [L]
t = time since the release of the solute [T]
Co = initial solute concentration [M/L^3]
v = average linear groundwater velocity (= q / porosity) [L/T]
D = longitudinal dispersion coefficient [L^2/T]

The Ogata-Banks solution was derived from the following initial and boundary conditions:
Initial condition: C(x,0) = C_i
Upper BC: Dirichlet-type (or concentration-type BC); C(0,t) = Co
Lower BC: Semi-infinite
So our Case 1 is comparable to the Ogata-Banks solution. Comparing Ogata-Banks to case 2 was just our curiosity.

Ref: Ogata, A., and Banks, R. B. 1961. A solution of the differential equation of longitudinal dispersion in porous media. U.S. Geological Survey Professional Paper 411-A, A1-A9.

Jirka
Posts: 5268
Joined: Sat Mar 16, 2002 3:47 pm
Location: USA
Location: Riverside, CA

### Re: Boundary condition implementation in Hydrus-1D

I think that the answer is obvious. The numerical model does not consider semi-infinite profile (as the analytical solution does); the soil profile is always finite and the boundary conditions are specified at a given locations and this BC affect the concentration profile in its near vicinity. If you want to compare the two solutions, you need to apply the BC (e.g., a zero gradient) at a location far from locations for which you compare the results. J.

Hue
Posts: 6
Joined: Tue Jul 16, 2019 12:32 pm
Location: Germany

### Re: Boundary condition implementation in Hydrus-1D

Thanks Jirka. That really solved the puzzle.