open TELEMAC-MASCARET The mathematically superior suite of solvers

Trying to get atmospheric exchange working. program crashed with diff3d non convergence and temperatures set to nan.

My model runs with no problems in 3D with tidal water level boundary and open flow boundaries.

I have tried to run with only surface heat exchange (Solar heating has been ignored for now). ATMOSEXCH=1 runs OK. With ATMOSEXCH=2 the run crashes in the diff3d with non convergence and nan's generated. I have looked at the values of ATABOS and BTABOS variables, and both ATMOSEXCH options give comparable values across the whole model. A constant water temperature is set (20C), and a constant air temperature is set (21C) in both cases. Solar heating of the water column is ignored at this stage.

Note the crash occurs after a number of timesteps. Up to that point temperatures look good. Changes to temperature are less than 0.1 degree. Expecting to see a slow increase in temperature.

I have used the example code in bord3d, and checked the code against a number of references.

I have been looking this over for a number of days with no luck.

Model options activated are DELWAQ output and therefore MASS LUMPING on. I have tried various diffusion options as well. In bord3d I forced (with ATMOSEXCH=2) BTABOS to zero - and model ran.

Running with telemac3d v7p1 under linux.

Suggestion - the exchange_with_atmosphere source code should specify the units of all required variables that may be input as timeseries.

Am I missing something important?

Comments appreciated. John

Hello,

After all your tests it looks like a numerical instability caused by BTABOS (explicit part of the boundary condition), so we would have to look at the formula with this option 2. If it is really a matter of stability, reducing the time-step should solve the problem or allow to go further. I am not sure that a flux condition at a boundary like nut*dC/dn=ATABOS*(C-C0)=ATABOS*C+BTABOS is always stable (I anticipate that it is this kind of condition that we have). It is probably unstable if C is too far from C0.

With best regards,

Jean-Michel Hervouet

What exactly do we need to define for BTABOS. Is it the end value of the variable value or the change to it. The current historical code assumes that it is the change, that is delT/s applied to that cell/element. This relates to my discussion below on documentation.

The input from heat exchange as programmed is in units of delTm/s. The code doesnt seem to do the conversion from energy/m^2 to delT per cell/element. Applying this correction (deviding by top cell depth/thickness) now leads to instant instability.

Looking through the code base trying to get an idea of the schematisation used in the diffusion code,and taking another source term TRAIN as an example, I am confused.
TRAIN source is added to the convected result (after converting input mass to concentration in the top element) in murd3d AND also added in the the diffusion calculation? This does not seem right.

There appears to be large changes in diff3d for mass lumping. Is something been missed for the surface boundary conditions treatment?

Is there a document describing the schematisation and discretisation used in the model for the various computational units? Many questions I have can be answered by a decent document. It is not very clear in the code whilst decyphering the vector manipulations. It seems that much of the geometry components are incorporated in the calculations at a lower level and not for example in the diff3d source.

Cheers
John

Hello,

(after a while... back from Cuba...)

Actually it is really like I said in my previous post :

nut*dC/dn=ATABOS*(C-C0)=ATABOS*C+BTABOS

where nut is the turbulent viscosity. The term nut*dC/dn is then integrated along the boundary to give a flux of matter leaving or entering the domain. The unit of nut is M**2/T (square of a distance divided by time). If I am not wrong (jetlag) the unit of ATABOS is then M/T, a velocity. BTABOS is then a velocity multipled by the unit of C.

With best regards,

Jean-Michel Hervouet

Sorry took a while to get back.

Thanks for that explanation. I took a slightly different appraoch and treated the surface exchange as a source rather than a boundary, which seems a logical approach.

Part of my problem was that the code documents these boundary conditions as dirichlet BCs yet arent they really neumann BCs. The manual does describe what you have said in the section on the boundary conditions file. I didnt link that to surface BC - sorry.

On the documentation side ...

Is there any (opensource) documentation that describes the numerical construct of the code eg including what the different integration routines do eg op=MASVEC? The manual references a Telemac3d technical note - I have not found this anywhere.

I will have another look at this is the near future.

A small suggestion. There are a number of routines that are called for quite different variables eg cvdf3d (and diff3d). It would be very useful to pass as an argument a tag that indicates which variables it is being operated on - I have done this in my copy, using a character string tag. It allows easier debugging
and also specific switches to be set if required.

My learning curve is improving. It is easy to get a base HD model working, with simple tracers etc. Adding in other 'important' features is a bit more involved (at the coding level) and debugging especially if running multi core is harder. Single core debug is hard because the program slows down so much. However the system IS comprehensive.

Thankyou
John

Hello,

For what concerns basic subroutines like MASVEC, did you look at the guide for programming in the Telemac system, which is available on this website ? In the section home/manuals..

opentelemac.org/index.php/component/jdow...mming-v6p0?Itemid=56

With best regards,

Jean-Michel Hervouet