open TELEMAC-MASCARET The mathematically superior suite of solvers

Dear Users,

I need your advice if you do not mind. I have a strongly stratified environment in my domain. All three boundaries are closed wall and one boundary is open. At open boundary I have various prescribed tracers values for each layer, which is different then in domain itself. The total number of horizontal layers is 30.

My question is: could I use one section of closed wall boundary for prescribing different (then in domain and at open boundary) tracers values. For example I wish it to represent a glacier front with temperature minus 1C and Salinity 0%, but not having an input of water from that boundary (still keep boundary as closed wall).

Sorry if I did not explain it properly. I wonder if it is possible and if you could please put me into right direction. I think I might need to make some changes in bord3d.f. Please advise. Thanks in advance.

Kind Regards!

Violeta

Hello,

A tracer cannot be imposed at a solid boundary, this would spoil mass conservation for example. What you can do is setting flux conditions across the solid boundaries, in the form (for a tracer T): dT/dn=coef*(T-T0), with coef sufficiently strong so that T becomes rapidly equal to T0, which means that behind your solid boundary there is a value T0 influencing T. T0 can vary on the vertical. On Telemac side you have viscosity*dT/dn=AFBOR*T+BFBOR, so you can find AFBOR and BFBOR giving you the right formula. AFBOR and BFOR can be given in the boundary conditions file, but in this case they are given per point and constant on the vertical, so in your case it would be better to give them in bord3d.f.

Anyway I let you think of it before going further, maybe there are other solutions,

With best regards,

Jean-Michel Hervouet

Dear Jean-Michel,

Thank you very much for your explanations. At the moment I am running case with modifying ATBOR and BTBOR in boudary file as constant on vertical profile. I will see the results. But you are right, ideally I will need it to be varied vertically at specific section at closed boundary.

I would appreciate if you could help me with code on bord3d subroutine. I already have section there coded for tracers vertical profile at open boundary. I probably will need to specify there the nodes number and plane number in fortran on which section i wish to make the change. If you could just help me with start and then I will experiment until I get it right, I still have to learn a lot in Fortran:). I have attached fortran file. Sorry. Thanks a mil. Only if you have time.

Kind Regards!

Violeta

Dear Jean-Michele/Users,

I have problem with specifying ATBOR and BTBOR in bord3d.f subroutine. I need to specify tracers profile behind closed wall at specific nodes only varying vertically. I have looked at complim.f and probably have to use similar to below script, maybe for each NPLAN, as i have them in total as 30. Could please hint that I am in right direction, or maybe give me more suggestions so I could experiment, as so far I am getting errors only. Thank you in advance.
The nodes I wish to change temperature and profile are 1265;1239;1222;1221;1219;1255;1304.

IF(TRAC.AND.IELMT.EQ.1) THEN
!
DO K=1,NPTFR
!
KP1=MESH%KP1BOR%I(K)
!
IF(KP1.NE.K) THEN
IF(LITBOR(K).EQ.LITBOR(KP1)) THEN
LITBOR(K+NPTFR) = LITBOR(K)
ELSEIF( LITBOR(K ).EQ.KLOG .OR.
& LITBOR(KP1).EQ.KLOG ) THEN
LITBOR(K+NPTFR) = KLOG
ELSEIF( LITBOR(K ).EQ.KADH .OR.
& LITBOR(KP1).EQ.KADH ) THEN
LITBOR(K+NPTFR) = KADH
ELSEIF( LITBOR(K ).EQ.KENTU .OR.
& LITBOR(KP1).EQ.KENTU ) THEN
LITBOR(K+NPTFR) = KENTU
ELSEIF( LITBOR(K ).EQ.KSORT .OR.
& LITBOR(KP1).EQ.KSORT ) THEN
LITBOR(K+NPTFR) = KSORT
ELSEIF( LITBOR(K ).EQ.KINC .OR.
& LITBOR(KP1).EQ.KINC ) THEN
LITBOR(K+NPTFR) = KINC
ELSE
WRITE(LU,*) 'CONDITION INITIALE QUADRATIQUE DE U ','K ',K,
& ' NON PREVUE POUR LIUBOR = ',LIUBOR(K),
& ' ET LIUBOR(K+1) = ',LIUBOR(MESH%KP1BOR%I(K))
CALL PLANTE(1)
STOP
ENDIF
TBOR(K+NPTFR) = (TBOR(K)+TBOR(KP1)) *0.5D0
ATBOR(K+NPTFR) = (ATBOR(K)+ATBOR(KP1))*0.5D0
BTBOR(K+NPTFR) = (BTBOR(K)+BTBOR(KP1))*0.5D0
ENDIF
!
ENDDO
ENDIF

Kind Regards!

Violeta

Hello Violetta,

Using KP1BOR only allows you to navigate on the 2D border on the horizontal. In your case,the coefficients to be changed are those dedicated to flux across lateral boundaries, for tracers, which are:

ATABOL and BTABOL

The nodes you want to change are given in global numbers, like 1265, so you will have to find their 2D boundary point number, let's say IPTFR, now assume you want to change the coefficients of the point on plane IPLAN above 2D boundary point IPTFR, for tracer ITRAC, the address is :

ATABOL%ADR(ITRAC)%P%R((IPLAN-1)*NPTFR2+IPTFR), where NPTFR2 is the number of boundary points in 2D (you go up one plane by adding NPTFR2). Same for BTABOL.

Once it works with hardcoded boundary points numbers, we can add a loop to retrieve these numbers from the global numbers, even in parallel.

I hope this is clear enough, though complicated, if not just ask again...

With best regards,

Jean-Michel Hervouet

Thank you very much Jean-Michele for suggestions and explanations. I have started working on it, so far with no success, but I have to spend more time on it and try harder. I will let you know tomorrow if I still have problems. Thanks again for your help.

Kind Regards!

Violeta

Dear Jean-Michele,

Thanks again for your directions. I managed to modify code below which does not give me any errors (I hope it is logical). I have two quick questions:

1) Is there a way to write the global points in one line as DO IPTFR = 1265, 1239, ... or I have to do similar to below (provided it is correct) for each point?

2) Assuming I have to do it for some NPLANE to vary vertically, then I need to specify various values for ATBOR and BTBOR coefficients?

Please let me know what you think on below and if it is correct then I will try to do Loop for each tracer for some levels. Sorry, I did not mentioned it but I am running in parallel, so hopefully it can be done in parallel too. Thanks a mil for your time.

IF(NTRAC.EQ.1) THEN
DO IPTFR = 1265,NPTFR2
IPOIN2 = NBOR2%I(IPTFR)
DO IPLAN = 1,NPLAN
IBORD = (IPLAN-1)*NPTFR2 + IPTFR
ATBOR(K+NPTFR2)=ATABOL%ADR(ITRAC)%P%R(IPLAN-1)*NPTFR2+IPTFR
BTBOR(K+NPTFR2)=BTABOL%ADR(ITRAC)%P%R(IPLAN-1)*NPTFR2+IPTFR
ENDDO
ENDDO
ENDIF
ELSEIF(NTRAC.EQ.2) THEN
DO IPTFR = 1265,NPTFR2
IPOIN2 = NBOR2%I(IPTFR)
DO IPLAN = 1,NPLAN
IBORD = (IPLAN-1)*NPTFR2 + IPTFR
ATBOR(K+NPTFR2)=ATABOL%ADR(ITRAC)%P%R(IPLAN-1)*NPTFR2+IPTFR
BTBOR(K+NPTFR2)=BTABOL%ADR(ITRAC)%P%R(IPLAN-1)*NPTFR2+IPTFR
ENDDO
ENDDO
ENDIF


Kind Regards!

Violeta

Hello,

I was not very clear, your values for the flux conditions across the solid walls must be put into ATABOL%ADR(ITRAC)%P%R((IPLAN-1)*NPTFR2+IPTFR) and BTABOL with the same address, so it is :

ATABOL%ADR(ITRAC)%P%R((IPLAN-1)*NPTFR2+IPTFR)= something
BTABOL%ADR(ITRAC)%P%R((IPLAN-1)*NPTFR2+IPTFR)= another something

that you need in your loop.

Then 1265 must be a boundary point number. If those numbers are in a sequence you can consider doing a loop, if not you can put their numbers in an array :

INTEGER ADR(10)
ADR(1)=1265
ADR(2)=...

and then:
DO II=1,10
IPTFR=ADR(II)
here what was in the DO IPTFR loop
ENDDO

AFBOR and BFBOR in my #14311 was just a generic name (given in diffusion subroutines) for such coefficients and in your case of tracers and lateral boundaries it is ATABOL and BTABOL.

With best regards,

JMH

Dear Jean-Michel,

Thanks again for your help and advice. I will do so.

Kind Regards!

Violeta

Dear Jean-Michel,

Sorry for bothering you again. Good news is that I did not get any error messages after compiling the case with fortran (in parallel on my cluster I have to do it separately before I submit it on the queue). However, after starting simulation, it crashes with error messages. Please see attached.
Could it be the issue with parallelism? Also i have attached Fortran file where I have loop with ATABOL etc. (it starts on line 607, as file is too big). Please advise what else I can try. Sorry if I made any mistakes again in writing the Loop, as that could be the cause of the problem. Thank you.

Kind Regards!

Violeta