open TELEMAC-MASCARET The mathematically superior suite of solvers

Hello,

In order to calculate history effects in time- and space-varying bed roughness (using subroutine ride_VR.f), I need to retain information at the old time step (specifically total roughness KS and ripple roughness KSR). This can be done quite easily in v7p2 using Private Arrays in Sisyphe, at least when a single processor is involved. However this approach seems to break down when multi-processors are used, apparently due to the complicated sequencing of the computations over the partitioned grid. This is certainly the case when 3-way coupling of T2D-Tom-Sis is attempted.

In T2D (v7p2, declarations_telemac2d.f) information from the old time step is obtainable, e.g.:
! U AND V AT OLD TIME-STEP
TYPE(BIEF_OBJ), TARGET :: UN,VN
with similar definitions: HN, TN, AKN, EPN at the old time step for: H, T, AK and EP, respectively. Also increments exist in T2D, e.g.
! INCREMENT OF H IN ONE TIME STEP
TYPE(BIEF_OBJ), TARGET :: DH
with similar increments: DU, DV and DHN for: U, V and HN, respectively.

However, there seem to be no old-time step or increment variables in declarations_sisyphe.f

So my question is whether anyone can advise how best to retain old-time step information in Sisyphe for multi-processing purposes?
Can old-time step variables, like those available in T2D, be safely introduced into Sisyphe by a reasonably competent user!?

Any advice or experiences would be appreciated.

Thanks, Alan Davies

Hello Alan,

If you have something that works in scalar mode with private user arrays I see no reason why it would not work in parallel, so this is perhaps the first point to investigate.

Otherwise you can create a KSN and a KSRN by declaring them:

TYPE(BIEF_OBJ) :: KSN,KSRN

BUT before being used they must be allocated by a CALL BIEF_ALLVEC(...) as is done for KS and KSR in point_sisyphe.f AND an important point is that this allocation must be done ONCE only (so it is best placed in an initialisation part of the code).

A SAVE statement at the beginning of the subroutine where they are declared will prevent forgetting them between two calls.

Then you can use them, e.g.:

CALL OS('X=Y ',X=KSN,Y=KS)


With best regards,

Jean-Michel

Hello Jean-Michel,

Thank you for your helpful reply. Your first comment suggests that I am writing to the private array in the wrong place and/or in the wrong way for multiprocessing. So I will look at this again, before creating new variables along the lines suggested (which may prove the best option).

What I have done hitherto, and successfully with 1 processor, is, at the end of RIDE_VR.f, to export an array of modified dune roughness values KSCD_NEW, I=1,NPOIN, to Unit 23, which is then rewound. The modification is a simple dune 'history effect' (which is described in a paper which Peter Robins and I have had accepted recently by Geomorphology). On returning to RIDE_VR.f at the next model time step these KSCD values are read back in from Unit 23 and are redesignated as KS_NOW, i.e. values of dune roughness to be used at the old time step in the history calculation. Unit 23 is rewound, and the process then repeats (calculation of new equilibrium dune roughness, comparison with KS_NOW, and application of the 'history effect'). The procedure has proved to be very robust.

However, where multiprocessing is concerned, I have noticed that the model subdomains are of unequal size (different NPOIN value for each subdomain) and that they are processed by RIDE_VR in a different order at successive model time steps. What then seems to go wrong with my procedure above is that information exported at the end of each call to RIDE_VR is somehow getting scrambled upon being read in at the start of the next call (possibly due to the way that I am exporting subdomain information, rewinding 23 and then reading from 23, but for the wrong sub-domain?!).

My 'incomplete' understanding of the multiprocessing sequence will by now fully apparent to you, and so I will write no more for now! But thank you very much for your reply which, hopefully, will steer me in the right direction.

With best wishes, Alan

Hello Jean-Michel,

As a closure on this, I am now using Private Arrays correctly, which was not the case before, and things are working well without any major changes to the code, as you anticipated!

Thanks again, Alan

Hello Jean-Michel,

Is it possible to define a BIEF_OBJ locally in a subroutine? I'm try to migrate the sediment sliding subroutine from SISYPHE to TELEMAC-3D because the high sedimentation on tidal flats.

I need to define a BIEF_OBJ so it could use PARCOM in parallel as it does in SISYPHE.

The other alternative way is to use PRIVE2D. But PRIVE2D contains 4 components. Can I use PARCOM for only one component in PRIVE2D?

Kind regards,
Qilong

Hi Qilong

I've already defined a BIEF_OBJ locally for different purpose so it works.

regards

Thanks for the help!

I will do it as you mentioned.

Kind regards,
Qilong

Hi Christophe,

I added the following part in FONVAS.f:

But I always got the errors saying

t3dfort.f(5510): error #6351: The number of subscripts is incorrect.
I1=MESH2D%IKLE%I(IELEM,1)

---

^
t3dfort.f(5511): error #6351: The number of subscripts is incorrect.
I2=MESH2D%IKLE%I(IELEM,2)

---

^
t3dfort.f(5512): error #6351: The number of subscripts is incorrect.
I3=MESH2D%IKLE%I(IELEM,3)

---

^
t3dfort.f(5514): error #6351: The number of subscripts is incorrect. [R]
X2=MESH2D%XEL%R(IELEM,2)

---

^
t3dfort.f(5515): error #6351: The number of subscripts is incorrect. [R]
X3=MESH2D%XEL%R(IELEM,3)

---

^
t3dfort.f(5516): error #6351: The number of subscripts is incorrect. [R]
Y2=MESH2D%YEL%R(IELEM,2)

---

^
t3dfort.f(5517): error #6351: The number of subscripts is incorrect. [R]
Y3=MESH2D%YEL%R(IELEM,3)

---

^

I'm wondering how to use the members IKLE, XEL and YEL in MESH2D?

Kind regards,
Qilong

Hi
You should check in the code but I think here IKLE is a vector of size ndp*nelem and not an array with 2 components ...
For XEL and YEL, as they are BIEF_OBJ, I will try to remove the %R.

Hope this helps

Thanks, Christophe!

In the Telemac developer guide, it says
But I checked the code. Indeed IKLE originally is a vector instead of a NELEM-by-3 matrix in the code. However, in some places, somehow it is reconstructed as a NELEM-by-3 matrix.

If it's a vector, is the order in IKLE like

[ELEM_1(1) ELEM_1(2) ELEM_1(3) ELEM_2(1) ELEM_2(2) ELEM_2(3) ...] ?

Kind regards,
Qilong