Welcome, Guest
Username: Password: Remember me

TOPIC: Initialize Salinity in 3D Mesh

Initialize Salinity in 3D Mesh 12 years 8 months ago #4112

  • qilong
  • qilong's Avatar
  • OFFLINE
  • Expert Boarder
  • Posts: 340
  • Thank you received: 33
Hello,

Last week I found my 3D hydrodynamics model became unstable when Salinity was added. Now I think I might found the reason. In my FORTRAN file, I use subroutine "CONDIM" to initialize salinity. Here is the code I wrote:

IF(NTRAC.GT.0) THEN

FILESAL='D:\LinuxExchange\Big_Mesh\Salinity.txt'
OPEN (LUSAL,FILE=FILESAL,FORM='FORMATTED',STATUS='OLD')

DO I=1,55855
READ (LUSAL,*) SA_NPOIN(I), SA_X(I), SA_Y(I), SA_S(I)
ENDDO

DO ITRAC=1,NTRAC
! CALL OS( 'X=C ', X=TA%ADR(ITRAC)%P, C=TRAC0(ITRAC))
IF (ITRAC .EQ. 1) THEN
DO I=1, NPOIN2
TA%ADR(ITRAC)%P%R(I) = SA_S(I)
WRITE(LU,*) I,SA_NPOIN(I),TA%ADR(ITRAC)%P%R(I)
ENDDO
ENDIF
ENDDO

ENDIF

I found if I disable the line "CALL OS( 'X=C ', X=TA%ADR(ITRAC)%P, C=TRAC0(ITRAC))", the model will become unstable. Then I looked at the result and found that the salinity was only imposed in the bottom layer and the other layers have initial values I gave in the steering file.

My question is how to impose salinity for all the layers? Which variable could be used in such case?

Thanks in advance!
The administrator has disabled public write access.

Re: Initialize Salinity in 3D Mesh 12 years 8 months ago #4113

  • jmhervouet
  • jmhervouet's Avatar
Hello,

The CALL OS correctly initializes all points at the value TRAC0(ITRAC), but your loop is done on NPOIN2, which is the number of points in the 2D mesh, you should do it on NPOIN2*NPLAN or MESH3D%NPOIN to reach all the layers. The CALL OS will do the loop on TA%ADR(ITRAC)%P%DIM1 which is the total number of 3D points.

With best regards,

Jean-Michel Hervouet
The administrator has disabled public write access.

Re: Initialize Salinity in 3D Mesh 12 years 8 months ago #4114

  • qilong
  • qilong's Avatar
  • OFFLINE
  • Expert Boarder
  • Posts: 340
  • Thank you received: 33
Thanks you very much!
It works now!
The administrator has disabled public write access.

Re: Initialize Salinity in 3D Mesh 11 years 10 months ago #6824

  • robinson
  • robinson's Avatar
Hello,

I have similar problem to initialize salinity in 3D mesh. having failed to initialize salinity with previous computation result automatically, I decide to do that manually with 2D salinity field and assign to each layer. But very strange is that my 'condim.f' works well in scalar mode, but when it comes to parallel mode, errors always come out(error.PNG). Here is the code:
      REWIND 27     
      DO J=1,NPOIN2
         READ(27,*) NPO(J),SA(J)
      ENDDO
!      
      DO ITRAC=1,NTRAC
!        CALL OS( 'X=C     ', X=TA%ADR(I)%P, C=TRAC0(I))
        IF(ITRAC .EQ. 1)THEN
           DO IPLAN=1,NPLAN
            DO J=1,NPOIN2
              IF(NCSIZE .GT. 0)THEN  
                  I=MESH3D%KNOGL%I((IPLAN-1)*NPOIN2+J)
                  IF(I.NE.0)THEN
                    TA%ADR(ITRAC)%P%R(I)=SA(J)
                  ENDIF
              ELSE
                 TA%ADR(ITRAC)%P%R((IPLAN-1)*NPOIN2+J)=SA(J)
              ENDIF
            ENDDO
          ENDDO
         ENDIF
       ENDDO
could you please give me a hint where the mistake is ..
thank you in advance!

with best regard.
Robinson
Attachments:
The administrator has disabled public write access.

Initialize Salinity in 3D Mesh 11 years 10 months ago #6825

  • pham
  • pham's Avatar
  • OFFLINE
  • Administrator
  • Posts: 1559
  • Thank you received: 602
Hello Robinson,

I think that there are at least 2 problems with your implementation.
The 1st one is that the FORMATTED DATA FILE 2 (corresponding to canal number 27) is not partitionned like a SERAFIN file in parallel (the last argument in telemac3dv6p2.dico is SCAL). Thus, the numbering of the fields is not changed even when running in parallel and it does not correspond to the general fields in TELEMAC (like velocities, tracers...). In particular, you have to read all the information for every node for this file, not only NPOIN2 (the number of 2D nodes of the mesh of one subdomain in parallel), e.g. NPOIN27 if it is the number of lines in your FORMATTED DATA FILE 2.

The 2nd one is related to the numbering of the nodes in parallel in your loops and the definition of MESh3D%KNOGL%I if this does exist (I have to check).

I would first try something like:
      REWIND 27     
      DO J=1,NPOIN27 ! OR HARD CODE OF THE NUMBER OF LINES IN FORMATTED DATA FILE 2
        READ(27,*) NPO(J),SA(J)
      ENDDO
!      
!      DO ITRAC=1,NTRAC
!        CALL OS( 'X=C     ', X=TA%ADR(I)%P, C=TRAC0(I))
!        IF(ITRAC.EQ.1) THEN
      ITRAC = 1
      DO IPLAN=1,NPLAN
        DO J=1,NPOIN2
          IF(NCSIZE.GT.0) THEN  
            TA%ADR(ITRAC)%P%R((IPLAN-1)*NPOIN2+J)=SA(MESH2D%KNOLG%I(J))
          ELSE
            TA%ADR(ITRAC)%P%R((IPLAN-1)*NPOIN2+J)=SA(J)
          ENDIF
        ENDDO
      ENDDO
!        ENDIF
!      ENDDO

Hope this helps,

Chi-Tuan
The administrator has disabled public write access.
The following user(s) said Thank You: robinson
Moderators: pham

The open TELEMAC-MASCARET template for Joomla!2.5, the HTML 4 version.