open TELEMAC-MASCARET The mathematically superior suite of solvers

Hy community,
I'm writing for a problem I can't solve by myself.
in my research, I need to simulate constructed wetland for the study of solute transport processes and contaminant removal efficiency. I have a very simple domain: an ideal wetland, rectangular shape 50mx200m, 0.5m water depht, inlet and outlet 10m wide. For having proper nominal residence time for the solute (6-7days), the entering discharge is very low, about 0.010m^3/s (and so, very low velocities in the domain).
My aim is to simulate the hydrodynamics, to reach the steady state (and mass conservation Qin=Qout), and then to use the steady-state flow condition for studying the solute fate. And do this for several bathymetric configurations.

The main problem in my simulation is that I can't reach the steady state!
I don't understand which is the problem...maybe the very low velocities?
Now I'm using
TIME STEP = 0.05
LISTING PRINTOUT PERIOD = 20
DURATION = 800000
STOP IF A STEADY STATE IS REACHED = YES
STOP CRITERIA = 1.E-5;1.E-5;1.E-4
and the simulation stops for reaching the time limit (800000s!)
but really, looking at the command prompt, it seems the Q values are very similar in between two printout periods (1s):

ITERATION 15100200 TIME: 8 D 17 H 43 MN 30.0001 S ( 755010.0001 S)

---

ADVECTION STEP

---

DIFFUSION-PROPAGATION STEP
EQUNOR (BIEF) : 0 ITERATIONS, ABSOLUTE PRECISION: 0.8930390E-04

---

BALANCE OF WATER VOLUME
VOLUME IN THE DOMAIN : 6861.298 M3
FLUX BOUNDARY 1: -0.1999181E-01 M3/S ( >0 : ENTERING <0 : EXITING )
FLUX BOUNDARY 2: 0.2000000E-01 M3/S ( >0 : ENTERING <0 : EXITING )
RELATIVE ERROR IN VOLUME AT T = 0.7550E+06 S : 0.2617471E-07
MAXIMUM COURANT NUMBER: 0.5767724E-02
TIME-STEP : 0.5000000E-01

================================================================================
ITERATION 15100220 TIME: 8 D 17 H 43 MN 31.0001 S ( 755011.0001 S)

---

ADVECTION STEP

---

DIFFUSION-PROPAGATION STEP
EQUNOR (BIEF) : 0 ITERATIONS, ABSOLUTE PRECISION: 0.9663666E-04

---

BALANCE OF WATER VOLUME
VOLUME IN THE DOMAIN : 6861.298 M3
FLUX BOUNDARY 1: -0.1999181E-01 M3/S ( >0 : ENTERING <0 : EXITING )
FLUX BOUNDARY 2: 0.2000000E-01 M3/S ( >0 : ENTERING <0 : EXITING )
RELATIVE ERROR IN VOLUME AT T = 0.7550E+06 S : 0.2911310E-07
MAXIMUM COURANT NUMBER: 0.5767724E-02
TIME-STEP : 0.5000000E-01

And it seems that if I use bigger steady state criteria:
STOP CRITERIA = 1.E-4;1.E-4;1.E-4
it stops after few iterations with Qout very different from Qin!

Does anybody give me help?
In the attached files, I uploaded a simplified example with: flat bathymetry, bigger discharge (0.020m3/s, but it also doesn't work), no tracer.
I'm only interested to have flow calculation working!

Thanks, Eleonora

Hi,

considering that I didn't have reply to my initial post, I tried to check something more by myself.

Even if the simulation doesn't stop for reaching steady-state criteria, I checked if mass conservation is satisfied after 800000s duration by comparing Qin and Qout (discharge at the inlet and outlet).
This is what I had:
1) from command prompt
flux boundary 2: 0.02000000 m^3/s (this is the inlet, where I imposed the discharge=0.02 m^3/s)
flux boundary 1: 0.01999972 m^3/s (this is the outlet, where I imposed the water depth)

2) from BlueKeneue (creating a field vector uh+vh and integrating along inlet and oulet lines)
Qin= 0.02000000 m^3/s
Qout=0.0187634 m^3/s

So, which is the correct Qout value? 0.01999972 m^3/s (so I could be satisfied by the mass conservation) or 0.0187634 m^3/s (about 93.8% of Qout)?

Waiting for some help about this last question, I'm trying to run the simulation for longer duration (1200000s!)...

Thanks, Eleonora

Dear Eleonora,


I don't know how BK computes a flux across a line. TELEMAC-2D gives you the good fluxes.

About the keyword "STOP IF A STEADY STATE IS REACHED": in the prompt can you see this messages:

THE STEADY STATE HAS BEEN REACHED
and
STEADY STATE REACHED, TELEMAC-2D IS STOPPED

?

I'm looking the code. I think previous messages are printed but the stop no!


best regards,
J

PS: you can restart calculation from a results file to reduce computational time.

Hi Jacques,
thanks for your reply.
in the prompt I only saw, when telemac stops,the messages:
"END OF TIME LOOP
CALL OF P_EXIT IN ITS VOID VERSION
CORRECT END OF RUN"
I don't know if somewhere there is a "the steady state has been reached".

if I can trust more to the values of Q out (flux boundary) showed in the prompt than the one calculated in BK, maybe the steady state is reached, with the exiting discharge almost equal to the entering discharge. I think I could be satisfied ...

Thanks for your suggestion: I ran a very long simulation (flat bathymetry, 1.200.000s) and reached Qout=0.0199997 m^3/s ... (not sure if steady state is reached inside the domain, considering Telemac always stops for end of time). And I decided to use this result for running other simulation with different bathymetry... even if the velocity field is different because of the different bathymetry, this seems faster for reaching Qout similar to Qin than starting from still water and Qout=0 !

Best. Eleonora

Dear Eleonora,

Steady state criterion seems very strict (you can have a look in $HOMETEL/source/telemac2d/steady.f subroutine, and with a dimention by the way). There are probably (reflecting) waves (due to initialisation, zero velocity) in your domain, that is why you do not reach the "steady.f state".

Could you test some iterations from a restart file with theses options:?

Best regards,
J

I'll try it, thanks!
Eleonora

Hi,
I tried Jacques's suggestion and I didn't notice big differences: the simulation reached the "end of time loop" with a value Qout similar to what I had in my simulations. So "steady state criteria" weren't satisfied.

Really, yesterday I discovered one mistake I made: in my geometry file, I added a roughness mesh (because in my old plans I wanted to test non homogeneous roughness)...and the value was about 35. (at the beginning I thought to use Chezy)
but now in the .cas file I impose Manning:
LAW OF BOTTOM FRICTION = 4
FRICTION COEFFICIENT = 0.025

1) in this contrasting case, does Telemac use roughness information found in .slf file, or in .cas file?

2) yesterday I started a new simulation with roughness mesh = 0.025 in the .slf file, to avoid contrasting info between .slf and .cas, but Qout I see in the prompt command (still running) seems much more instable than in the uncorrected simulation!
Do you suggest to delete roughness mesh from the geometry file?
Do you suggest to use another law for bottom friction? or another value for Manning?

In attachment my files, again the simpler case with flat bathymetry and no tracer, only for calculating hydrodinamics. (starting from t=0, still water, because it is the first run I make with updated roughness).

Thanks, Eleonora

Dear Eleonora,

I did not want to be rude against steady.f. I think you should not use it.

(with a time step = 1e-14s, I suppose you reach "steady.f" in 1 iteration...)

You can follow evolution in time of 2 or 3 probes (far from boundary conditions) to check if steady state is reached.

About your questions:
1) .slf file
2) As long your roughness is constant in space, I suggest to remove from the .slf file for sake of simplicity. I am not surprised by your (unwanted) numerical experience, I think steady state is more difficult to obtain with a smooth bottom.

Best regards,
J

Hi Jacques,
I'll try your suggestions:
- add in the .cas
TYPE OF ADVECTION = 5;5
SCHEME OPTION FOR ADVECTION OF VELOCITIES = 4
- delete steady-state stopping criteria
- use longer simulation duration and check steady state using some probes.
--> can I add this in the .cas and see the values in the listing printout, or do I check them in BlueKenue from T2DRES file sometimes during the simulation?
- delete unused roughness field in the selafin

I need to start from a flat bottom in order to compare its results (hydrodynamics and solute transport) with other bathymetric configurations...

Really, It doesn't seem a complicated case to solve: rectangular domain, flat bottom, constant roughness, Qin constant. The only issue seems the long duration because very slow regime...
I'm considering to start with a completely new set up in the file .cas, maybe after so many trial tests I can't see that there is something wrong affecting my simulation!... (I wonder what...)

I'll try again. Thanks for your help!
Eleonora

Dear Eleonora,

Your test case is not so simple.
You have to carefully initialise the flow to reduce computational time.
1) start from a file (restart)
2) start with an initialised velocity and a slope for the inlet flowrate to reach Qin smoothly

Else waves are generated. These waves cross the domain (good! Thx to non diffusive numerical schemes - 5 is more diffusive than 1 but preserve mass -) but are reflected at the outlet (H is prescribed... - btw you can try Thompson bcs for the outlet: OPTION FOR LIQUID BOUNDARIES)

You can add probes in .cas (LIST OF POINTS and NAMES OF POINTS) or follow probes in BK or extract with the embedded python tools a point in your results file (run_selafin.py —help).

Best regards,
Jacques