OGS5 Water and Tracer Pulse Unsaturated Richards Flow

Hey all,

I’m new to OGS and work with OGS 5 and coupled BRNS for my master thesis. For starters I try to implement a simple example of a tracer and simultaneous water pulse into an unsaturated 2D domain in OGS 5.

My Set-Up (in short):

The domain represents a 0.5m x 0.5m cross section through a soil column of a sandy soil, characterised by van Genuchten parameters.

IC is a PRESSURE1 GRADIENT (cross section is a XZ plane) with -31600 Pa at the bottom boundary and a gradient of 9810 Pa/m upwards (thus -36505 Pa at the upper boundary) and a conc. of 0 mmol C/g soil tracer in the entire domain.

BC for Richards flow PRESSURE1 are a constant bottom boundary of -31600 Pa and a time variable source term entering from the upper boundary. The source term represents a two hour long water pulse simulated by a Gauss curve in the .rfd file which adds a total of 10E-05 m3 of water into the domain. Additionally, a tracer is added throughout the first hour (also by a Gauss curve in .rfd) by a time dependant upper boundary condition for the tracer.

My Problem:

By running this set-up I get considerably large negative concentrations of tracer (in the order of E-02) one to two nodes from the concentration peak.

My question concerns the set-up of the numeric solver. Unfortunately, I couldn’t figure out which solver setting is appropriate for this problem. I tried to apply the numeric and time set-up of the HM benchmark in \DENSITY-DEPENDENT_FLOW\goswami-clement\wholeBC_PressAsHead_tri as it as well combines MASS_TRANSPORT with RICHARDS_FLOW. However, from the available documentation I couldn’t figure out how this set-up works and if it is applicable for my problem. Could anyone please help me in this regard - or might there be any other mistakes in my OGS set-up? I pasted the .num and .tim file below (not allowed to uplaod anything as a new user). Please let me know if you required any other files. I would highly appreciate any help and/or feedback!

Cheers,

Erik


NUM

$OVERALL_COUPLING
   2 12
#NUMERICS
  $PCS_TYPE
   RICHARDS_FLOW
  $LINEAR_SOLVER			; linear solver (2=BiCGSTab) , convergence crit.
   2 6 1e-14 1000 1 100 4
  $NON_LINEAR_ITERATION		; solver type , convergence crit. , iteration no. , preconditioner , error tollerance
   PICARD ERNORM 25 0 1e-10
  $COUPLING_CONTROL
   ERNORM 0.001
#NUMERICS
  $PCS_TYPE
   MASS_TRANSPORT
  $LINEAR_SOLVER
   2 6 1e-14 5000 1 100 4
  $COUPLING_CONTROL
   ERNORM 0.001
#STOP

TIM

#TIME_STEPPING
  $PCS_TYPE
   RICHARDS_FLOW
  $TIME_START
   0.0
  $TIME_END
   200
  $TIME_UNIT
   DAYS
  $TIME_CONTROL
   SELF_ADAPTIVE
   3 1.3
   7 0.8
   MIN_TIME_STEP
   1e-05
   MAX_TIME_STEP
   10
   INITIAL_STEP_SIZE			; set initial time step below Pulse temporal resolution
   0.001
   ITERATIVE_TYPE
   COUPLED
   STAY
   3
#TIME_STEPPING
  $PCS_TYPE
   MASS_TRANSPORT
  $TIME_START
   0.0
  $TIME_END
   200
  $TIME_UNIT
   DAYS
  $TIME_CONTROL
   SELF_ADAPTIVE
   3 1.3
   7 0.8
   MIN_TIME_STEP
   1e-05
   MAX_TIME_STEP
   10
   INITIAL_STEP_SIZE			; set initial time step below Pulse temporal resolution
   0.001
   ITERATIVE_TYPE
   COUPLED
   STAY
   3
#STOP

Current output at time 4E-02 (during pulse application). Node lenght is 0.01 in X and Z.

1 Like

Hi Erik,

sounds like a numerical oscillation at the concentration front. Your mesh seems quite coarse, try to refine it, especially at the top.
Can you provide the remaining input files, please? Maybe as a download link to DropBox or similar sevice.

Best,
Johannes

Hey Johannes,

thanks for the hint! I also experienced oscillating behaviours in other set-ups.
Is there a recommended spatial resolution/ rule of thumb?

Files are uploaded to DropBox now:

I refined the temporal resolution to make sure that time step size with the SELF_ADAPTIVE setting does not exceed temporal discretisation of the input pulses while they are active.

Cheers,
Erik

The SELF_ADAPTIVE time control scheme respects the CFL number (velocitytime_step/element_size!<=1) to ensure stability. However, your element size is too large. I refined the mesh around the top using GMSH and got rid of the oscillation. You’ll have to play around with it to find the optimal balance between accuracy and speed. You’ll find the file to define the mesh in GMSH under ~/mesh/.geo. You gotta use the ogs-data-explorer to convert GMSH *.msh files to ogs *.msh. For reactive transport simulations you may want to list the final time steps used in the conservative tracer transport simulation from the log.file to use them for a fixed time stepping scheme to save computation time.

Files are here:

1 Like

Thank you very much Johannes!
Works perfectly, also with my other simulations where I had problems with oscillating behaviours. I will familiarise a bit more with GMSH and try to optimise computational speed. Also thanks for the hint with fixed time stepping there!

Cheers, Erik

Short addition: OGS Data Explorer can’t handle *.msh generated by most recent GMSH release (version 4.3.0). Don’t know how far this goes back, but it works for me with GMSH version 3.0.6 downloadable from http://gmsh.info/bin/ .

Revisiting this topic again, using

$MASS_LUMPING
1

substantally reduced oscillatory behaviour. Are there any major drawbacks with this approach or is it an legimate alternative to higher spatial resolution?