w3triamd Module Reference

Reads triangle and unstructured grid information. More...

Functions/Subroutines
subroutine	readmsh (NDS, FNAME)
	Reads triangle and unstructured grid information from GMSH files. More...
subroutine	readmsh_iobp (NDS, FNAME)
	Reads triangle and unstructured grid information from GMSH files. More...
subroutine	get_boundary_status (STATUS)
	Boundary status (code duplication). More...
subroutine	readmshobc (NDS, FNAME, TMPSTA, UGOBCOK)
	Reads open boundary information for UNST grids following GMESH type format. More...
subroutine	ug_getopenboundary (TMPSTA, ZBIN, ZLIM)
	Defines open boundary points based on depth. More...
subroutine	spatial_grid
	Calculates triangle areas and reorders the triangles to have them oriented counterclockwise. More...
subroutine	nvectri
	Calculate cell tools: inward normal, angles and length of edges. More...
subroutine	count (TRIGPTEMP)
	Calculate global and maximum number of connection for array allocations. More...
subroutine	coordmax
	Calculate first point and last point coordinates, and minimum and maximum edge length. More...
subroutine	area_si (IMOD)
	Define optimized connection arrays (points and triangles) for spatial propagation schemes. More...
subroutine	is_in_ungrid (IMOD, XTIN, YTIN, ITOUT, IS, JS, RW)
	Determine whether a point is inside or outside an unstructured grid, and returns index of triangle and interpolation weights. More...
subroutine	is_in_ungrid2 (IMOD, XTIN, YTIN, FORCE, ITOUT, IS, JS, RW)
	Determine whether a point is inside or outside an unstructured grid, and returns index of triangle and interpolation weights. More...
subroutine	ug_gradients (PARAM, DIFFX, DIFFY)
	Calculate gradients at a point via its connection. More...
subroutine	w3nestug (DISTMIN, FLOK)
	UGTYPE nesting initialization. More...
subroutine	set_iobp (MASK, STATUS)
	Setup boundary pointer. More...
subroutine	get_boundary (MNP, MNE, TRIGP, IOBP, NEIGHBOR_PREV, NEIGHBOR_NEXT)
	Find boundary points. More...
subroutine	triang_indexes (I, INEXT, IPREV)
	Set indices of the triangle. More...
subroutine	set_ug_iobp ()
	Redefines the values of the boundary points and angle pointers based on the MAPSTA array. More...
subroutine	fix_periodcity (I1, I2, I3, XGRD, YGRD, PT)
	Adjust element longitude coordinates for elements straddling the dateline with distance of ~360 degrees. More...

Variables
integer, save	n_outside_boundary
integer, dimension(:), allocatable, save	outside_boundary

Detailed Description

Reads triangle and unstructured grid information.

Look for namelist with name NAME in unit NDS and read if found.

Author

F. Ardhuin

A. Roland

Date

26-Jan-2014

Function/Subroutine Documentation

area_si()

subroutine w3triamd::area_si

(

integer, intent(in)

IMOD

)

Define optimized connection arrays (points and triangles) for spatial propagation schemes.

The storage is optimize especially considering the iterative solver used. The schemes used are vertex-centered, a point has to be considered within its median dual cell. For a given point, the surface of the dual cell is one third of the sum of the surface of connected triangles.

This routine is from WWM developed in Darmstadt(Aaron Roland).

Parameters

[in]

IMOD

Model number to point to.

Author

A. Roland

Date

23-Aug-2011

Definition at line 1337 of file w3triamd.F90.

    !/ -------------------------------------------------------------------
    !/                  +-----------------------------------+
    !/                  | WAVEWATCH III           NOAA/NCEP |
    !/                  |           A. Roland               |
    !/                  |                        FORTRAN 90 |
    !/                  | Last update :          23-Aug-2011|
    !/                  +-----------------------------------+
    !/
    !/    15-May-2007 : Origination: adjustment from the WWM code       ( version 3.13 )
    !/    23-Aug-2011 : Removes double entries in VNEIGH                ( version 4.04 )
    !/
    !
    !  1. Purpose :
    !
    !      Define optimized connection arrays (points and triangles) for spatial propagation schemes.
    !
    !  2. Method :
    !
    !  3. Parameters :
    !
    !  4. Subroutines used :
    !
    !  5. Called by :
    !
    !       Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !      READTRI    Subr. Internal  Unstructured mesh definition.
    !     ----------------------------------------------------------------
    !
    !  6. Error messages :
    !
    !  7. Remarks :
    !
    !     The storage is optimize especially considering the iterative solver used.
    !     The schemes used are vertex-centered, a point has to be considered within its
    !     median dual cell. For a given point, the surface of the dual cell is one third
    !     of the sum of the surface of connected triangles.
    !     This routine is from WWM developped in Darmstadt(Aaron Roland)
    !
    !  8. Structure :
    !
    !  9. Switches :
    !
    ! 10. Source code :
    !
    !/ ------------------------------------------------------------------- /
 
    USE w3gdatmd
#ifdef W3_S
    USE w3servmd, ONLY: strace
#endif
    IMPLICIT NONE
    !/ input
 
    INTEGER, INTENT(IN) :: IMOD
 
    !/ local parameters
 
    INTEGER :: COUNTER,ifound,alreadyfound
    INTEGER :: I, J, K, II
    INTEGER :: IP, IE, POS, POS_I, POS_J, POS_K, IP_I, IP_J, IP_K
    INTEGER :: I1, I2, I3, IP2, CHILF(NX)
    INTEGER :: TMP(NX), CELLVERTEX(NX,COUNTRI,2)
    INTEGER :: COUNT_MAX
    DOUBLE PRECISION   :: TRIA03
    INTEGER, ALLOCATABLE :: PTABLE(:,:)
 
#ifdef W3_S
    INTEGER                      ::  IENT = 0
#endif
    !/ ------------------------------------------------------------------- /
 
#ifdef W3_S
    CALL strace (ient, 'AREA_SI')
#endif
 
    si(:) = 0.d0
    !
    ccon(:) = 0     ! Number of connected Elements
    DO ie = 1, ntri
      i1 = trigp(1,ie)
      i2 = trigp(2,ie)
      i3 = trigp(3,ie)
      ccon(i1) = ccon(i1) + 1
      ccon(i2) = ccon(i2) + 1
      ccon(i3) = ccon(i3) + 1
      tria03 = 1./3. * tria(ie)
      si(i1) = si(i1) + tria03
      si(i2) = si(i2) + tria03
      si(i3) = si(i3) + tria03
    ENDDO
 
    cellvertex(:,:,:) = 0 ! Stores for each node the Elementnumbers of the connected Elements
    ! and the Position of the Node in the Element Index
 
    chilf = 0
 
    DO ie = 1, ntri
      DO j=1,3
        i = trigp(j,ie)!INE(J,IE)
        chilf(i) = chilf(i)+1
        cellvertex(i,chilf(i),1) = ie
        cellvertex(i,chilf(i),2) = j
      END DO
    ENDDO
    !
    ! Second step in storage, the initial 3D array CELLVERTEX, is transformed in a 1D array
    ! the global index is J . From now, all the computation step based on these arrays must
    ! abide by the conservation of the 2 loop algorithm (points + connected triangles)
    !
    index_cell(1)=1
    j = 0
    DO ip = 1, nx
      DO i = 1, ccon(ip)
        j = j + 1
        ie_cell(j)  = cellvertex(ip,i,1)
        pos_cell(j) = cellvertex(ip,i,2)
      END DO
      index_cell(ip+1)=j+1
    END DO
 
    IF (.NOT. fsnimp) RETURN
 
    j = 0
    DO ip = 1, nx
      DO i = 1, ccon(ip)
        j = j + 1
      END DO
    END DO
 
    count_max = j
 
    ALLOCATE(ptable(count_max,7))
 
    j = 0
    ptable(:,:) = 0.
    DO ip = 1, nx
      DO i = 1, ccon(ip)
        j = j + 1
        ie    = ie_cell(j)
        pos   = pos_cell(j)
        i1 = trigp(1,ie)
        i2 = trigp(2,ie)
        i3 = trigp(3,ie)
        IF (pos == 1) THEN
          pos_j = 2
          pos_k = 3
        ELSE IF (pos == 2) THEN
          pos_j = 3
          pos_k = 1
        ELSE
          pos_j = 1
          pos_k = 2
        END IF
        ip_i = ip
        ip_j = trigp(pos_j,ie)
        ip_k = trigp(pos_k,ie)
        ptable(j,1) = ip_i
        ptable(j,2) = ip_j
        ptable(j,3) = ip_k
        ptable(j,4) = pos
        ptable(j,5) = pos_j
        ptable(j,6) = pos_k
        ptable(j,7) = ie
      END DO
    END DO
 
    !        WRITE(*,'("+TRACE......",A)') 'SET UP SPARSE MATRIX POINTER ... COUNT NONZERO ENTRY'
 
    j = 0
    k = 0
    DO ip = 1, nx
      tmp(:) = 0
      DO i = 1, ccon(ip)
        j = j + 1
        ip_j  = ptable(j,2)
        ip_k  = ptable(j,3)
        pos   = ptable(j,4)
        tmp(ip)   = 1
        tmp(ip_j) = 1
        tmp(ip_k) = 1
      END DO
      k = k + sum(tmp)
    END DO
 
    nnz => grids(imod)%NNZ
    nnz = k
 
    !          WRITE(*,'("+TRACE......",A)') 'SET UP SPARSE MATRIX POINTER ... SETUP POINTER'
 
    ALLOCATE (grids(imod)%JAA(nnz))
    ALLOCATE (grids(imod)%IAA(nx+1))
    ALLOCATE (grids(imod)%POSI(3,count_max))
    jaa   => grids(imod)%JAA
    iaa   => grids(imod)%IAA
    posi  => grids(imod)%POSI
 
    j = 0
    k = 0
    iaa(1) = 1
    jaa    = 0
    DO ip = 1, nx ! Run through all rows
      tmp = 0
      DO i = 1, ccon(ip)         ! Check how many entries there are ...
        j = j + 1                ! this is the same J index as in IE_CELL
        ip_j  = ptable(j,2)
        ip_k  = ptable(j,3)
        tmp(ip)   = 1
        tmp(ip_j) = 1
        tmp(ip_k) = 1
      END DO
      DO i = 1, nx               ! Run through all columns
        IF (tmp(i) .GT. 0) THEN  ! this is true only for the connected points
          k = k + 1
          jaa(k) = i
        END IF
      END DO
      iaa(ip + 1) = k + 1
    END DO
 
    posi = 0
    j = 0
    DO ip = 1, nx
      DO i = 1, ccon(ip)
        j = j + 1
        ip_j  = ptable(j,2)
        ip_k  = ptable(j,3)
        DO k = iaa(ip), iaa(ip+1) - 1
          IF (ip   == jaa(k)) posi(1,j)  = k
          IF (ip_j == jaa(k)) posi(2,j)  = k
          IF (ip_k == jaa(k)) posi(3,j)  = k
          IF (k == 0) THEN
            WRITE(*,*) .EQ.'ERROR IN AREA_SI K  0'
            stop
          END IF
        END DO
      END DO
    END DO
 
    DEALLOCATE(ptable)
 

References w3gdatmd::ccon, w3gdatmd::fsnimp, w3gdatmd::grids, w3gdatmd::iaa, w3gdatmd::ie_cell, w3gdatmd::index_cell, w3gdatmd::jaa, w3gdatmd::nnz, w3gdatmd::ntri, w3gdatmd::pos_cell, w3gdatmd::posi, w3gdatmd::si, w3servmd::strace(), w3gdatmd::tria, and w3gdatmd::trigp.

Referenced by readmsh(), and w3initmd::w3init().

coordmax()

subroutine w3triamd::coordmax

Calculate first point and last point coordinates, and minimum and maximum edge length.

Author

F. Ardhuin

Date

15-May-2008

Definition at line 1249 of file w3triamd.F90.

    !/ -------------------------------------------------------------------
    !/                  +-----------------------------------+
    !/                  | WAVEWATCH III           NOAA/NCEP |
    !/                  |           F. Ardhuin              |
    !/                  |                        FORTRAN 90 |
    !/                  | Last update :          15-May-2008|
    !/                  +-----------------------------------+
    !/
    !/    15-May-2007 : Origination.                        ( version 3.13 )
    !/
    !  1. Purpose :
    !
    !      Calculate first point and last point coordinates, and minimum and maximum edge length.
    !
    !  2. Method :
    !
    !  3. Parameters :
    !
    !  4. Subroutines used :
    !
    !  5. Called by :
    !
    !       Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !      READTRI    Subr. Internal  Unstructured mesh definition.
    !     ----------------------------------------------------------------
    !
    !  6. Error messages :
    !
    !  7. Remarks :
    !
    !  8. Structure :
    !
    !  9. Switches :
    !
    ! 10. Source code :
    !
    !/ ------------------------------------------------------------------- /
    USE w3gdatmd
#ifdef W3_S
    USE w3servmd, ONLY: strace
#endif
    IMPLICIT NONE
#ifdef W3_S
    INTEGER                      ::  IENT = 0
#endif
 
 
#ifdef W3_S
    CALL strace (ient, 'COORDMAX')
#endif
    !
    ! maximum of coordinates s
    !
    maxx = maxval(xgrd(1,:))
    maxy = maxval(ygrd(1,:))
    !
    ! minimum of coordinates
    !
    x0 = minval(xgrd(1,:))
    y0 = minval(ygrd(1,:))
    !
    !maximum and minimum length of edges
    !
    dxymax = maxval(len(:,:))
    sx = minval(len(:,:))
    sy = sx
    !

References w3gdatmd::dxymax, w3gdatmd::len, w3gdatmd::maxx, w3gdatmd::maxy, w3servmd::strace(), w3gdatmd::sx, w3gdatmd::sy, w3gdatmd::x0, w3gdatmd::xgrd, w3gdatmd::y0, and w3gdatmd::ygrd.

Referenced by readmsh(), and w3initmd::w3init().

count()

subroutine w3triamd::count

(

integer, dimension(:,:), intent(in)

TRIGPTEMP

)

Calculate global and maximum number of connection for array allocations.

Parameters

[in]

TRIGPTEMP

Temporary array of triangle vertices.

Author

A. Roland

F. Ardhuin

Date

15-May-2008

Definition at line 1135 of file w3triamd.F90.

    !/ -------------------------------------------------------------------
    !/                  +-----------------------------------+
    !/                  | WAVEWATCH III           NOAA/NCEP |
    !/                  |           A. Roland               |
    !/                  |           F. Ardhuin              |
    !/                  |                        FORTRAN 90 |
    !/                  | Last update :          15-May-2008|
    !/                  +-----------------------------------+
    !/
    !/    15-May-2007 : Origination.                        ( version 3.13 )
    !/
    !
    !  1. Purpose :
    !
    !      Calculate global and maximum number of connection for array allocations .
    !
    !  2. Method :
    !
    !  3. Parameters :
    !     Parameter list
    !     ----------------------------------------------------------------
    !       NTRI         Int.   I   Total number of triangle.
    !       TRIGPTEMP    Int    I   Temporary array of triangle vertices
    !       COUNTRI      Int    O   Maximum number of connected triangle
    !                               for a given points
    !       COUNTOT      Int    O   Global number of triangle connection
    !                               for the whole grid.
    !     ----------------------------------------------------------------
    !  4. Subroutines used :
    !
    !  5. Called by :
    !
    !       Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !      READTRI    Subr. Internal  Unstructured mesh definition.
    !     ----------------------------------------------------------------
    !
    !  6. Error messages :
    !
    !  7. Remarks :
    !
    !  8. Structure :
    !
    !  9. Switches :
    !
    ! 10. Source code :
    !
    !/ ------------------------------------------------------------------- /
    USE w3gdatmd
#ifdef W3_S
    USE w3servmd, ONLY: strace
#endif
    IMPLICIT NONE
 
 
    !/ parameter list
 
    INTEGER,INTENT(IN) :: TRIGPTEMP(:,:)
    !/ ------------------------------------------------------------------- /
    !/ local parameter
 
    INTEGER               :: CONN(NX)
    INTEGER               :: COUNTER, IP, IE, I, J, N(3)
#ifdef W3_S
    INTEGER                      ::  IENT = 0
#endif
    !/------------------------------------------------------------------------
 
#ifdef W3_S
    CALL strace (ient, 'COUNT')
#endif
 
    countri=0
    countot=0
    conn(:)= 0
 
    !
    !calculate the number of connected triangles for a given point.
    !
 
    DO ie = 1,ntri
      n(:) = 0.
      n(1) = trigptemp(1,ie)
      n(2) = trigptemp(2,ie)
      n(3) = trigptemp(3,ie)
      conn(n(1)) = conn(n(1)) + 1
      conn(n(2)) = conn(n(2)) + 1
      conn(n(3)) = conn(n(3)) + 1
    ENDDO
 
    countri = maxval(conn)
    !
    ! calculate the global number of connections available through the mesh
    !
    j=0
    DO  ip=1,nx
      DO i=1,conn(ip)
        j=j+1
      ENDDO
    ENDDO
    countot=j
 

References w3gdatmd::countot, w3gdatmd::countri, w3gdatmd::ntri, and w3servmd::strace().

Referenced by readmsh().

fix_periodcity()

subroutine w3triamd::fix_periodcity	(	integer, intent(in)	I1,
		integer, intent(in)	I2,
		integer, intent(in)	I3,
		double precision, dimension(:,:), intent(in)	XGRD,
		double precision, dimension(:,:), intent(in)	YGRD,
		real*8, dimension(3,2), intent(out)	PT
	)

Adjust element longitude coordinates for elements straddling the dateline with distance of ~360 degrees.

Detect if element has nodes on both sides of dateline and adjust coordinates so that all nodes have the same sign.

Parameters

[in]	I1
[in]	I2
[in]	I3
[in]	XGRD
[in]	YGRD
[out]	PT

Author

Steven Brus

Ali Abdolali

Date

21-May-2020

Definition at line 3059 of file w3triamd.F90.

    !/
    !/                  +-----------------------------------+
    !/                  | WAVEWATCH III           NOAA/NCEP |
    !/                  |        Steven Brus                |
    !/                  |        Ali Abdolali               |
    !/                  |                        FORTRAN 90 |
    !/                  | Last update :         21-May-2020 |
    !/                  +-----------------------------------+
    !/
    !/    21-May-2020 : Origination.                        ( version 6.07 )
    !/
    !/
    !  1. Purpose :
    !
    !     Adjust element longitude coordinates for elements straddling the
    !     dateline with distance of ~360 degrees
    !
    !  2. Method :
    !
    !     Detect if element has nodes on both sides of dateline and adjust
    !     coordinates so that all nodes have the same sign
    !
    !  3. Parameters :
    !
    !     Parameter list
    !     ----------------------------------------------------------------
    IMPLICIT NONE
    INTEGER, INTENT(IN) :: I1, I2, I3
    DOUBLE PRECISION, INTENT(IN) :: XGRD(:,:), YGRD(:,:)
    real*8, INTENT(OUT) :: pt(3,2)
    !     ----------------------------------------------------------------
    !
    !     Local variables.
    !     ----------------------------------------------------------------
    INTEGER :: I
    INTEGER :: R1GT180, R2GT180, R3GT180
    !     ----------------------------------------------------------------
    !
    !  4. Subroutines used :
    !
 
    !  5. Called by :
    !
    !      Name         Type  Module   Description
    !     ----------------------------------------------------------------
    !      SPATIAL_GRID Subr. W3TRIAM  Triangle area calculation
    !      NVECTRI      Subr. W3TRIAM  Edge length, angle, normal calcuation
    !      IS_IN_UNGRID Subr. W3TRIAM  Point in element calculation
    !     ----------------------------------------------------------------
    !
    !  6. Error messages :
    !
    !       None.
    !
    !  7. Remarks :
    !
    !  8. Structure :
    !
    !  9. Switches :
    !
    ! 10. Source code :
    !/ ------------------------------------------------------------------- /
 
    pt(1,1) = xgrd(1,i1)
    pt(1,2) = ygrd(1,i1)
    pt(2,1) = xgrd(1,i2)
    pt(2,2) = ygrd(1,i2)
    pt(3,1) = xgrd(1,i3)
    pt(3,2) = ygrd(1,i3)
 
 
    r1gt180 = merge(1, 0, abs(pt(3,1)-pt(2,1)).GT.180)
    r2gt180 = merge(1, 0, abs(pt(1,1)-pt(3,1)).GT.180)
    r3gt180 = merge(1, 0, abs(pt(2,1)-pt(1,1)).GT.180)
    ! if R1GT180+R2GT180+R3GT180 .eq. 0 the element does not cross the dateline
    ! if R1GT180+R2GT180+R3GT180 .eq. 1 the element contains the pole
    ! if R1GT180+R2GT180+R3GT180 .eq. 2 the element crosses the dateline
 
 
    IF ( r1gt180 + r2gt180 == 2 ) THEN
      pt(3,1)=pt(3,1)-sign(360.0d0,(pt(3,1)-pt(2,1)))
    ELSE IF ( r2gt180 + r3gt180 == 2 ) THEN
      pt(1,1)=pt(1,1)-sign(360.0d0,(pt(1,1)-pt(2,1)))
    ELSE IF ( r1gt180 + r3gt180 == 2 ) THEN
      pt(2,1)=pt(2,1)-sign(360.0d0,(pt(2,1)-pt(3,1)))
    ENDIF
 
    RETURN

Referenced by is_in_ungrid(), nvectri(), and spatial_grid().

get_boundary()

subroutine w3triamd::get_boundary	(	integer, intent(in)	MNP,
		integer, intent(in)	MNE,
		integer, dimension(3,mne), intent(in)	TRIGP,
		integer*2, dimension(mnp), intent(inout)	IOBP,
		integer, dimension(mnp), intent(inout)	NEIGHBOR_PREV,
		integer, dimension(mnp), intent(inout)	NEIGHBOR_NEXT
	)

Find boundary points.

Parameters

[in]	MNP
[in]	MNE
[in]	TRIGP
[in,out]	IOBP
[in,out]	NEIGHBOR_PREV
[in,out]	NEIGHBOR_NEXT

Author

Aron Roland

Mathieu Dutour-Sikiric

Date

01-Jun-2018

Definition at line 2481 of file w3triamd.F90.

    !/
    !/                  +-----------------------------------+
    !/                  | WAVEWATCH III           NOAA/NCEP |
    !/                  |                                   |
    !/                  | Aron Roland (BGS IT&E GmbH)       |
    !/                  | Mathieu Dutour-Sikiric (IRB)      |
    !/                  |                                   |
    !/                  |                        FORTRAN 90 |
    !/                  | Last update :        01-June-2018 |
    !/                  +-----------------------------------+
    !/
    !/    01-June-2018 : Origination.                        ( version 6.04 )
    !/
    !  1. Purpose : find boundary points
    !  2. Method :
    !  3. Parameters :
    !
    !     Parameter list
    !     ----------------------------------------------------------------
    !     ----------------------------------------------------------------
    !
    !  4. Subroutines used :
    !
    !      Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !      STRACE    Subr. W3SERVMD Subroutine tracing.
    !     ----------------------------------------------------------------
    !
    !  5. Called by :
    !
    !      Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !     ----------------------------------------------------------------
    !
    !  6. Error messages :
    !  7. Remarks
    !  8. Structure :
    !  9. Switches :
    !
    !     !/S  Enable subroutine tracing.
    !
    ! 10. Source code :
    !
    !/ ------------------------------------------------------------------- /
#ifdef W3_S
    USE w3servmd, ONLY: strace
#endif
    !
    USE w3servmd, ONLY: extcde
    IMPLICIT NONE
    !/
    !/ ------------------------------------------------------------------- /
    !/ Parameter list
    !/
    !/ ------------------------------------------------------------------- /
    !/ Local PARAMETERs
    !/
#ifdef W3_S
    INTEGER, SAVE           :: IENT = 0
#endif
 
    INTEGER, INTENT(IN)             :: MNP, MNE, TRIGP(3,MNE)
    INTEGER*2, INTENT(INOUT)        :: IOBP(MNP)
    INTEGER, INTENT(INOUT)          :: NEIGHBOR_PREV(MNP)
    INTEGER, INTENT(INOUT)          :: NEIGHBOR_NEXT(MNP)
 
    INTEGER, POINTER :: STATUS(:)
    INTEGER, POINTER :: COLLECTED(:)
    INTEGER, POINTER :: NEXTVERT(:)
    INTEGER, POINTER :: PREVVERT(:)
 
    INTEGER :: IE, I, IP, IP2, IP3
    INTEGER :: ISFINISHED, INEXT, IPREV, ISTAT
    INTEGER :: IPNEXT, IPPREV, ZNEXT, ZPREV
#ifdef W3_S
    CALL strace (ient, 'GET_BOUNDARY')
#endif
    ALLOCATE(status(mnp), stat=istat)
    check_alloc_status( istat )
    ALLOCATE(collected(mnp), stat=istat)
    check_alloc_status( istat )
    ALLOCATE(prevvert(mnp), stat=istat)
    check_alloc_status( istat )
    ALLOCATE(nextvert(mnp), stat=istat)
    check_alloc_status( istat )
    neighbor_next = 0
    neighbor_prev = 0
    !  Now computing the next items
    status = 0
    nextvert = 0
    prevvert = 0
 
    DO ie=1,mne
      DO i=1,3
        CALL triang_indexes(i, inext, iprev)
        ip=trigp(i,ie)
        ipnext=trigp(inext,ie)
        ipprev=trigp(iprev,ie)
        IF (status(ip).EQ.0) THEN
          status(ip)=1
          prevvert(ip)=ipprev
          nextvert(ip)=ipnext
        END IF
      END DO
    END DO
    status(:)=0
    DO
      collected(:)=0
      DO ie=1,mne
        DO i=1,3
          CALL triang_indexes(i, inext, iprev)
          ip=trigp(i,ie)
          ipnext=trigp(inext,ie)
          ipprev=trigp(iprev,ie)
          IF (status(ip).EQ.0) THEN
            znext=nextvert(ip)
            IF (znext.EQ.ipprev) THEN
              collected(ip)=1
              nextvert(ip)=ipnext
              IF (nextvert(ip).EQ.prevvert(ip)) THEN
                status(ip)=1
              END IF
            END IF
          END IF
        END DO
      END DO
 
      isfinished=1
      DO ip=1,mnp
        IF ((collected(ip).EQ.0).AND.(status(ip).EQ.0)) THEN
          status(ip)=-1
          neighbor_next(ip)=nextvert(ip)
        END IF
        IF (status(ip).EQ.0) THEN
          isfinished=0
        END IF
      END DO
      IF (isfinished.EQ.1) THEN
        EXIT
      END IF
    END DO
 
    !  Now computing the prev items
    status = 0
    nextvert = 0
    prevvert = 0
    DO ie=1,mne
      DO i=1,3
        CALL triang_indexes(i, inext, iprev)
        ip=trigp(i,ie)
        ipnext=trigp(inext,ie)
        ipprev=trigp(iprev,ie)
        IF (status(ip).EQ.0) THEN
          status(ip)=1
          prevvert(ip)=ipprev
          nextvert(ip)=ipnext
        END IF
      END DO
    END DO
    status(:)=0
    DO
      collected(:)=0
      DO ie=1,mne
        DO i=1,3
          CALL triang_indexes(i, inext, iprev)
          ip=trigp(i,ie)
          ipnext=trigp(inext,ie)
          ipprev=trigp(iprev,ie)
          IF (status(ip).EQ.0) THEN
            zprev=prevvert(ip)
            IF (zprev.EQ.ipnext) THEN
              collected(ip)=1
              prevvert(ip)=ipprev
              IF (prevvert(ip).EQ.nextvert(ip)) THEN
                status(ip)=1
              END IF
            END IF
          END IF
        END DO
      END DO
 
      isfinished=1
      DO ip=1,mnp
        IF ((collected(ip).EQ.0).AND.(status(ip).EQ.0)) THEN
          status(ip)=-1
          neighbor_prev(ip)=prevvert(ip)     ! new code
        END IF
        IF (status(ip).EQ.0) THEN
          isfinished=0
        END IF
      END DO
      IF (isfinished.EQ.1) THEN
        EXIT
      END IF
    END DO
    !  Now making checks
    DO ip=1,mnp
      ip2=neighbor_next(ip)
      IF (ip2.GT.0) THEN
        ip3=neighbor_prev(ip2)
        IF (abs(ip3 - ip).GT.0) THEN
          WRITE(*,*) 'IP=', ip, ' IP2=', ip2, ' IP3=', ip3
          WRITE(*,*) 'We have a dramatic inconsistency'
          stop
        END IF
      END IF
    END DO
    !   Now assigning the boundary IOBP array
    DO ip=1,mnp
      IF (status(ip).EQ.-1 .AND. iobp(ip) .EQ. 1) THEN
        iobp(ip)=0
      END IF
    END DO
 
    DEALLOCATE(status, stat=istat)
    check_dealloc_status( istat )
    DEALLOCATE(collected, stat=istat)
    check_dealloc_status( istat )
    DEALLOCATE(nextvert, stat=istat)
    check_dealloc_status( istat )
    DEALLOCATE(prevvert, stat=istat)
    check_dealloc_status( istat )
 

References w3servmd::extcde(), w3servmd::strace(), and triang_indexes().

Referenced by set_iobp().

get_boundary_status()

subroutine w3triamd::get_boundary_status

(

integer*2, dimension(nx), intent(out)

STATUS

)

Boundary status (code duplication).

Parameters

[out]

STATUS

Author

Mathieu Dutour-Sikiric

Aron Roland

Date

01-May-2018

Definition at line 511 of file w3triamd.F90.

    !/
    !/                  +-----------------------------------+
    !/                  | WAVEWATCH III           NOAA/NCEP |
    !/                  |                                   |
    !/                  | Mathieu Dutour-Sikiric (IRB)      |
    !/                  | Aron Roland (BGS IT&E GmbH)       |
    !/                  |                                   |
    !/                  |                        FORTRAN 90 |
    !/                  | Last update :         01-Mai-2018 |
    !/                  +-----------------------------------+
    !/
    !/    01-Mai-2018 : Origination.                        ( version 6.04 )
    !/
    !  1. Purpose : boundary status (code duplication)
    !  2. Method :
    !  3. Parameters :
    !
    !     Parameter list
    !     ----------------------------------------------------------------
    !     ----------------------------------------------------------------
    !
    !  4. Subroutines used :
    !
    !      Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !      STRACE    Subr. W3SERVMD Subroutine tracing.
    !     ----------------------------------------------------------------
    !
    !  5. Called by :
    !
    !      Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !     ----------------------------------------------------------------
    !
    !  6. Error messages :
    !  7. Remarks
    !  8. Structure :
    !  9. Switches :
    !
    !     !/S  Enable subroutine tracing.
    !
    ! 10. Source code :
    !
    !/ ------------------------------------------------------------------- /
#ifdef W3_S
    USE w3servmd, ONLY: strace
#endif
    !
 
#ifdef W3_PDLIB
    use yowelementpool, only: ne_global
    use yownodepool, only: np_global
#endif
    USE w3gdatmd, ONLY : trigp, ntri, nx
    IMPLICIT NONE
    !/
    !/ ------------------------------------------------------------------- /
    !/ Parameter list
    !/
    !/ ------------------------------------------------------------------- /
    !/ Local PARAMETERs
    !/
#ifdef W3_S
    INTEGER, SAVE           :: IENT = 0
#endif
    !/
    !/ ------------------------------------------------------------------- /
    !/
    !
    integer*2, intent(out) :: STATUS(NX)
    INTEGER :: COLLECTED(NX), NEXTVERT(NX), PREVVERT(NX)
    INTEGER :: ISFINISHED, INEXT, IPREV
    INTEGER :: IPNEXT, IPPREV, ZNEXT, IP, I, IE
#ifdef W3_S
    CALL strace (ient, 'VA_SETUP_IOBPD')
#endif
    status(:) = 0
    DO ie=1,ntri
      DO i=1,3
        IF (i.EQ.1) THEN
          iprev=3
        ELSE
          iprev=i-1
        END IF
        IF (i.EQ.3) THEN
          inext=1
        ELSE
          inext=i+1
        END IF
        ip=trigp(i,ie)
        ipnext=trigp(inext,ie)
        ipprev=trigp(iprev,ie)
        IF (status(ip).EQ.0) THEN
          status(ip)=1
          prevvert(ip)=ipprev
          nextvert(ip)=ipnext
        END IF
      END DO
    END DO
    status(:)=0
    DO
      collected(:)=0
      DO ie=1,ntri
        DO i=1,3
          IF (i.EQ.1) THEN
            iprev=3
          ELSE
            iprev=i-1
          END IF
          IF (i.EQ.3) THEN
            inext=1
          ELSE
            inext=i+1
          END IF
          ip=trigp(i,ie)
          ipnext=trigp(inext,ie)
          ipprev=trigp(iprev,ie)
          IF (status(ip).eq.0) THEN
            znext=nextvert(ip)
            IF (znext.eq.ipprev) THEN
              collected(ip)=1
              nextvert(ip)=ipnext
              IF (nextvert(ip).eq.prevvert(ip)) THEN
                status(ip)=1
              END IF
            END IF
          END IF
        END DO
      END DO
      isfinished=1
      DO ip=1,nx
        IF ((collected(ip).eq.0).and.(status(ip).eq.0)) THEN
          status(ip)=-1
        END IF
        IF (status(ip).eq.0) THEN
          isfinished=0
        END IF
      END DO
      IF (isfinished.eq.1) THEN
        EXIT
      END IF
    END DO

References yowelementpool::ne_global, yownodepool::np_global, w3gdatmd::ntri, w3gdatmd::nx, w3servmd::strace(), and w3gdatmd::trigp.

is_in_ungrid()

subroutine w3triamd::is_in_ungrid	(	integer, intent(in)	IMOD,
		double precision, intent(in)	XTIN,
		double precision, intent(in)	YTIN,
		integer, intent(out)	ITOUT,
		integer, dimension(4), intent(out)	IS,
		integer, dimension(4), intent(out)	JS,
		real, dimension(4), intent(out)	RW
	)

Determine whether a point is inside or outside an unstructured grid, and returns index of triangle and interpolation weights.

This is the analogue for triangles of the FUNCTION W3GRMP.

Using barycentric coordinates defined as the ratio of triangle algebric areas which are positive or negative. Computes the 3 interpolation weights for each triangle until they are all positive.

Parameters

[in]	IMOD	Model number to point to.
[in]	XTIN	X-coordinate of target point.
[in]	YTIN	Y-coordinate of target point.
[out]	ITOUT	Model number to point to.
[out]	IS	I indices of vertices of enclosing grid cell.
[out]	JS	J indices of vertices of enclosing grid cell.
[out]	RW	Array of interpolation weights.

Author

Mathieu Dutour Sikiric, IRB

Aron Roland, Z&P

Fabrice Ardhuin

Date

26-Jan-2014

Definition at line 1605 of file w3triamd.F90.

    !/ -------------------------------------------------------------------
    !/                  +-----------------------------------+
    !/                  | WAVEWATCH III           NOAA/NCEP |
    !/                  |      Mathieu Dutour Sikiric, IRB  |
    !/                  |                 Aron Roland, Z&P  |
    !/                  |             Fabrice Ardhuin       |
    !/                  |                        FORTRAN 90 |
    !/                  | Last update :          26-Jan-2014|
    !/                  +-----------------------------------+
    !/
    !/ Adapted from other subroutine
    !/    15-Oct-2007 : Origination.                        ( version 3.13 )
    !/    21-Sep-2012 : Uses same interpolation as regular  ( version 4.08 )
    !/    26-Jan-2014 : Correcting bug in RW                ( version 4.18 )
    !/
    !  1. Purpose :
    !
    !      Determine whether a point is inside or outside an unstructured grid,
    !      and returns index of triangle and interpolation weights
    !      This is the analogue for triangles of the FUNCTION W3GRMP
    !
    !  2. Method :
    !
    !     Using barycentric coordinates defined as the ratio of triangle algebric areas
    !     which are positive or negative.
    !     Computes the 3 interpolation weights for each triangle until they are all positive
    !
    !  3. Parameters :
    !
    !     Parameter list
    !     ----------------------------------------------------------------
    !       IMOD    Int.   I   Model number to point to.
    !       XTIN    Real   I   X-coordinate of target point.
    !       YTIN    Real   I   Y-coordinate of target point.
    !       ITOUT    Int.   I   Model number to point to.
    !       IS,JS   I.A.   O   (I,J) indices of vertices of enclosing grid cell.
    !       RW      R.A.   O   Array of interpolation weights.
    !     ----------------------------------------------------------------
    !
    !  4. Subroutines used :
    !
    !     None
    !
    !  5. Called by :
    !
    !     WMGLOW, W3IOPP, WMIOPP, WW3_GINT
    !
    !  6. Error messages :
    !
    !     - Error checks on previous setting of variable.
    !
    !  7. Remarks :
    !
    !  8. Structure :
    !
    !  9. Switches :
    !
    !     !/S    Enable subroutine tracing.
    !     !/T    Enable test output
    !
    ! 10. Source code :
    !
    !  2. Method :
    !
    !     Using barycentric coordinates. Each coefficient depends on the mass of its related point in the interpolation.
    !
    !  3. Parameters :
    !
    !  4. Subroutines used :
    !
    !  5. Called by :
    !
    !       Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !      W3IOPP    Subr. Internal  Preprocessing of point output.
    !     ----------------------------------------------------------------
    !
    !  6. Error messages :
    !
    !  7. Remarks :
    !
    !      This subroutine is adjusted from CREST code (Fabrice Ardhuin)
    !      For a given output point, the algorithm enable to glance through all the triangles
    !      to find the one the point belong to, and then make interpolation.
    !
    !  8. Structure :
    !
    !  9. Switches :
    !
    !       !/LLG   Spherical grid.
    !       !/XYG   Carthesian grid.
    !
    ! 10. Source code :
    !
    !/ ------------------------------------------------------------------- /
    USE w3gdatmd
    USE w3servmd, ONLY: extcde
#ifdef W3_S
    USE w3servmd, ONLY: strace
#endif
    USE w3odatmd, ONLY: ndse
    IMPLICIT NONE
 
    !/ ------------------------------------------------------------------- /
    ! Parameter list
 
    INTEGER, INTENT(IN)            :: IMOD
    DOUBLE PRECISION, INTENT(IN)   :: XTIN, YTIN
    INTEGER, INTENT(OUT)           :: itout
    INTEGER, INTENT(OUT)           :: IS(4), JS(4)
    REAL, INTENT(OUT)              :: RW(4)
    !/ ------------------------------------------------------------------- /
    !local parameters
 
    DOUBLE PRECISION             :: x1, x2, x3
    DOUBLE PRECISION             :: y1, y2, y3
    DOUBLE PRECISION             :: s1, s2, s3, sg1, sg2, sg3
    real*8                       :: pt(3,2)
    INTEGER                      :: ITRI
    INTEGER                      :: I1, I2, I3
    INTEGER                      :: nbFound
#ifdef W3_S
    INTEGER                      ::  IENT = 0
    CALL strace (ient, 'IS_IN_UNGRID')
#endif
 
    !
    itout = 0
    nbfound=0
    itri = 0
    DO WHILE (nbfound.EQ.0.AND.itri.LT.grids(imod)%NTRI)
      itri = itri +1
      i1=grids(imod)%TRIGP(1,itri)
      i2=grids(imod)%TRIGP(2,itri)
      i3=grids(imod)%TRIGP(3,itri)
 
      CALL fix_periodcity(i1,i2,i3,grids(imod)%XGRD,grids(imod)%YGRD,pt)
      ! coordinates of the first vertex A
      x1 = pt(1,1)
      y1 = pt(1,2)
      ! coordinates of the 2nd vertex B
      x2 = pt(2,1)
      y2 = pt(2,2)
      !coordinates of the 3rd vertex C
      x3 = pt(3,1)
      y3 = pt(3,2)
      !with M = (XTIN,YTIN) the target point ...
      !vector product of AB and AC
      sg3=(y3-y1)*(x2-x1)-(x3-x1)*(y2-y1)
      !vector product of AB and AM
      s3=(ytin-y1)*(x2-x1)-(xtin-x1)*(y2-y1)
      !vector product of BC and BA
      sg1=(y1-y2)*(x3-x2)-(x1-x2)*(y3-y2)
      !vector product of BC and BM
      s1=(ytin-y2)*(x3-x2)-(xtin-x2)*(y3-y2)
      !vector product of CA and CB
      sg2=(y2-y3)*(x1-x3)-(x2-x3)*(y1-y3)
      !vector product of CA and CM
      s2=(ytin-y3)*(x1-x3)-(xtin-x3)*(y1-y3)
      IF ((s1*sg1.GE.0).AND.(s2*sg2.GE.0).AND.(s3*sg3.GE.0)) THEN
        itout=itri
        nbfound=nbfound+1
        is(1)=i1
        is(2)=i2
        is(3)=i3
        is(4)=1
        js(:)=1
        rw(1)=s1/sg1
        rw(2)=s2/sg2
        rw(3)=1.-rw(1)-rw(2)  !s3/sg3
        rw(4)=0.
      END IF
    ENDDO

References w3servmd::extcde(), fix_periodcity(), w3gdatmd::grids, w3odatmd::ndse, and w3servmd::strace().

Referenced by w3iopomd::w3iopp(), wmgridmd::wmglow(), and wmiopomd::wmiopp().

is_in_ungrid2()

subroutine w3triamd::is_in_ungrid2	(	integer, intent(in)	IMOD,
		double precision, intent(in)	XTIN,
		double precision, intent(in)	YTIN,
		integer, intent(in)	FORCE,
		integer, intent(out)	ITOUT,
		integer, dimension(4), intent(out)	IS,
		integer, dimension(4), intent(out)	JS,
		real, dimension(4), intent(out)	RW
	)

Determine whether a point is inside or outside an unstructured grid, and returns index of triangle and interpolation weights.

This is the analogue for triangles of the FUNCTION W3GRMP.

Using barycentric coordinates defined as the ratio of triangle algebric areas which are positive or negative. Computes the 3 interpolation weights for each triangle until they are all positive.

Parameters

[in]	IMOD	Model number to point to.
[in]	XTIN	X-coordinate of target point.
[in]	YTIN	Y-coordinate of target point.
[in]	FORCE
[out]	ITOUT	Model number to point to.
[out]	IS	I indices of vertices of enclosing grid cell.
[out]	JS	J indices of vertices of enclosing grid cell.
[out]	RW	Array of interpolation weights.

Author

Mathieu Dutour Sikiric, IRB

Aron Roland, Z&P

Fabrice Ardhuin

Date

26-Jan-2014

Definition at line 1807 of file w3triamd.F90.

    !/ -------------------------------------------------------------------
    !/                  +-----------------------------------+
    !/                  | WAVEWATCH III           NOAA/NCEP |
    !/                  |      Mathieu Dutour Sikiric, IRB  |
    !/                  |                 Aron Roland, Z&P  |
    !/                  |             Fabrice Ardhuin       |
    !/                  |                        FORTRAN 90 |
    !/                  | Last update :          26-Jan-2014|
    !/                  +-----------------------------------+
    !/
    !/ Adapted from other subroutine
    !/    15-Oct-2007 : Origination.                        ( version 3.13 )
    !/    21-Sep-2012 : Uses same interpolation as regular  ( version 4.08 )
    !/    26-Jan-2014 : Correcting bug in RW                ( version 4.18 )
    !/
    !  1. Purpose :
    !
    !      Determine whether a point is inside or outside an unstructured grid,
    !      and returns index of triangle and interpolation weights
    !      This is the analogue for triangles of the FUNCTION W3GRMP
    !
    !  2. Method :
    !
    !     Using barycentric coordinates defined as the ratio of triangle algebric areas
    !     which are positive or negative.
    !     Computes the 3 interpolation weights for each triangle until they are all positive
    !
    !  3. Parameters :
    !
    !     Parameter list
    !     ----------------------------------------------------------------
    !       IMOD    Int.   I   Model number to point to.
    !       XTIN    Real   I   X-coordinate of target point.
    !       YTIN    Real   I   Y-coordinate of target point.
    !       ITOUT    Int.   I   Model number to point to.
    !       IS,JS   I.A.   O   (I,J) indices of vertices of enclosing grid cell.
    !       RW      R.A.   O   Array of interpolation weights.
    !     ----------------------------------------------------------------
    !
    !  4. Subroutines used :
    !
    !     None
    !
    !  5. Called by :
    !
    !     WMGLOW, W3IOPP, WMIOPP, WW3_GINT
    !
    !  6. Error messages :
    !
    !     - Error checks on previous setting of variable.
    !
    !  7. Remarks :
    !
    !  8. Structure :
    !
    !  9. Switches :
    !
    !     !/S    Enable subroutine tracing.
    !     !/T    Enable test output
    !
    ! 10. Source code :
    !
 
 
    !  2. Method :
    !
    !     Using barycentric coordinates. Each coefficient depends on the mass of its related point in the interpolation.
    !
    !  3. Parameters :
    !
    !  4. Subroutines used :
    !
    !  5. Called by :
    !
    !       Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !      W3IOPP    Subr. Internal  Preprocessing of point output.
    !     ----------------------------------------------------------------
    !
    !  6. Error messages :
    !
    !  7. Remarks :
    !
    !      This subroutine is adjusted from CREST code (Fabrice Ardhuin)
    !      For a given output point, the algorithm enable to glance through all the triangles
    !      to find the one the point belong to, and then make interpolation.
    !
    !  8. Structure :
    !
    !  9. Switches :
    !
    !       !/LLG   Spherical grid.
    !       !/XYG   Carthesian grid.
    !
    ! 10. Source code :
    !
    !/ ------------------------------------------------------------------- /
    USE w3gdatmd
    USE w3servmd, ONLY: extcde
#ifdef W3_S
    USE w3servmd, ONLY: strace
#endif
    USE w3odatmd, ONLY: ndse
    IMPLICIT NONE
 
    !/ ------------------------------------------------------------------- /
    ! Parameter list
 
    INTEGER, INTENT(IN)            :: IMOD, FORCE
    DOUBLE PRECISION, INTENT(IN)            :: XTIN, YTIN
    INTEGER, INTENT(OUT)           :: itout
    INTEGER, INTENT(OUT)           :: IS(4), JS(4)
    REAL, INTENT(OUT)              :: RW(4)
    !/ ------------------------------------------------------------------- /
    !local parameters
 
    DOUBLE PRECISION             :: x1, x2, x3, D1, D2, D3, DISTMIN, DDMIN
    DOUBLE PRECISION             :: s1, s2, s3, sg1, sg2, sg3, smin, ssum
    DOUBLE PRECISION             :: y1, y2, y3
    INTEGER                      :: ITRI, ITRIS
    INTEGER                      :: I1, I2, I3
    INTEGER                      :: nbFound
    LOGICAL                      :: MAPSTAOK
#ifdef W3_S
    INTEGER                      ::  IENT = 0
    CALL strace (ient, 'IS_IN_UNGRID2')
#endif
 
    !
    itout = 0
    nbfound=0
    itri = 0
    itris = 1
    ssum = 0
    smin = 0
    DO WHILE (nbfound.EQ.0.AND.itri.LT.grids(imod)%NTRI)
      itri = itri +1
      i1=grids(imod)%TRIGP(1,itri)
      i2=grids(imod)%TRIGP(2,itri)
      i3=grids(imod)%TRIGP(3,itri)
      ! coordinates of the first vertex A
      x1=grids(imod)%XGRD(1,i1)
      y1=grids(imod)%YGRD(1,i1)
      ! coordinates of the 2nd vertex B
      x2=grids(imod)%XGRD(1,i2)
      y2=grids(imod)%XGRD(1,i2)
      !coordinates of the 3rd vertex C
      x3=grids(imod)%XGRD(1,i3)
      y3=grids(imod)%YGRD(1,i3)
      !with M = (XTIN,YTIN) the target point ...
      !vector product of AB and AC
      sg3=(y3-y1)*(x2-x1)-(x3-x1)*(y2-y1)
      !vector product of AB and AM
      s3=(ytin-y1)*(x2-x1)-(xtin-x1)*(y2-y1)
      !vector product of BC and BA
      sg1=(y1-y2)*(x3-x2)-(x1-x2)*(y3-y2)
      !vector product of BC and BM
      s1=(ytin-y2)*(x3-x2)-(xtin-x2)*(y3-y2)
      !vector product of CA and CB
      sg2=(y2-y3)*(x1-x3)-(x2-x3)*(y1-y3)
      !vector product of CA and CM
      s2=(ytin-y3)*(x1-x3)-(xtin-x3)*(y1-y3)
      !     ssum = ABS(s1*sg1)+ABS(s2*sg2)+ABS(s3*sg3)
      mapstaok = ((grids(imod)%MAPSTA(1,i1).GE.1).AND. &
           (grids(imod)%MAPSTA(1,i2).GE.1).AND.(grids(imod)%MAPSTA(1,i3).GE.1))
      IF (force.LT.2) mapstaok =.true.
      ssum = (xtin-(x1+x2+x3)/3.)**2+(ytin-(y1+y2+y2)/3.)**2
      IF (smin.EQ.0.AND. mapstaok ) smin=ssum
      !WRITE(6,*) 'ssum',ITRI,MAPSTAOK,ssum,smin
      IF (ssum.LT.smin .AND.  mapstaok  ) THEN
        smin=ssum
        itris=itri
      ENDIF
      IF ((s1*sg1.GE.0).AND.(s2*sg2.GE.0).AND.(s3*sg3.GE.0)) THEN
        itout=itri
        nbfound=nbfound+1
        is(1)=i1
        is(2)=i2
        is(3)=i3
        is(4)=1
        js(:)=1
        rw(1)=s1/sg1
        rw(2)=s2/sg2
        rw(3)=1.-rw(1)-rw(2)  !s3/sg3
        rw(4)=0.
      END IF
    ENDDO
    IF (itout.EQ.0.AND.force.GT.0) THEN
      itri=itris
      i1=grids(imod)%TRIGP(1,itri)
      i2=grids(imod)%TRIGP(2,itri)
      i3=grids(imod)%TRIGP(3,itri)
      ! coordinates of the first vertex A
      x1=grids(imod)%XGRD(1,i1)
      y1=grids(imod)%YGRD(1,i1)
      ! coordinates of the 2nd vertex B
      x2=grids(imod)%XGRD(1,i2)
      y2=grids(imod)%YGRD(1,i2)
      !coordinates of the 3rd vertex C
      x3=grids(imod)%XGRD(1,i3)
      y3=grids(imod)%YGRD(1,i3)
      d1=(xtin-x1)**2+(ytin-y1)**2
      d2=(xtin-x2)**2+(ytin-y2)**2
      d3=(xtin-x3)**2+(ytin-y3)**2
      IF (d1.LE.d2.AND.d1.LE.d3) is(1)=i1
      IF (d2.LE.d1.AND.d2.LE.d3) is(1)=i2
      IF (d3.LE.d2.AND.d3.LE.d1) is(1)=i3
      is(2:4)=1
      js(:)=1
      rw(1)=1
      rw(2:4)=0.
      itout=itri
 
    ENDIF

References w3servmd::extcde(), w3gdatmd::grids, w3odatmd::ndse, and w3servmd::strace().

nvectri()

subroutine w3triamd::nvectri

Calculate cell tools: inward normal, angles and length of edges.

To get inward pointing normals, triangle are glanced through anti-clockwisely.

Author

A. Roland

Date

15-May-2008

Definition at line 1004 of file w3triamd.F90.

    !/ -------------------------------------------------------------------
    !/                  +-----------------------------------+
    !/                  | WAVEWATCH III           NOAA/NCEP |
    !/                  |           A. Roland               |
    !/                  |                        FORTRAN 90 |
    !/                  | Last update :          15-May-2008|
    !/                  +-----------------------------------+
    !/
    !/    15-May-2007 : Origination: adjustment from the WWM code       ( version 3.13 )
    !/
    !
    !  1. Purpose :
    !
    !      Calculate cell tools: inward normal, angles and length of edges.
    !
    !  2. Method :
    !      To get inward pointing normals, triangle are glanced through anti-clockwisely
    !
    !
    !  3. Parameters :
    !
    !  4. Subroutines used :
    !
    !  5. Called by :
    !
    !       Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !      READTRI    Subr. Internal  Unstructured mesh definition.
    !     ----------------------------------------------------------------
    !
    !  6. Error messages :
    !
    !  7. Remarks :
    !
    !  8. Structure :
    !
    !  9. Switches :
    !
    ! 10. Source code :
    !
    !/ ------------------------------------------------------------------- /
    USE w3gdatmd
#ifdef W3_S
    USE w3servmd, ONLY: strace
#endif
    USE constants
 
    IMPLICIT NONE
    !
    !local parameter
    !
    INTEGER :: IP, IE
    INTEGER :: I1, I2, I3, I11, I22, I33
    !
    real*8    :: p1(2), p2(2), p3(2)
    real*8    :: r1(2), r2(2), r3(2)
    real*8    :: n1(2), n2(2), n3(2)
    real*8    :: tmp(3)
    real*8    :: tmpinv(3)
    real*8    :: pt(3,2)
#ifdef W3_S
    INTEGER                      ::  IENT = 0
#endif
    !/ ------------------------------------------------------------------- /
#ifdef W3_S
    CALL strace (ient, 'NVECTRI')
#endif
    !
    DO ie = 1, ntri
      !
      ! vertices
      !
      i1 = trigp(1,ie)
      i2 = trigp(2,ie)
      i3 = trigp(3,ie)
 
      CALL fix_periodcity(i1,i2,i3,xgrd,ygrd,pt)
 
      p1(1) = pt(1,1)
      p1(2) = pt(1,2)
      p2(1) = pt(2,1)
      p2(2) = pt(2,2)
      p3(1) = pt(3,1)
      p3(2) = pt(3,2)
      !
      ! I1 -> I2, I2 -> I3, I3 -> I1 (anticlockwise orientation is preserved)
      !
      r1 = p3-p2
      r2 = p1-p3
      r3 = p2-p1
 
      n1(1) = (-r1(2))
      n1(2) = ( r1(1))
      n2(1) = (-r2(2))
      n2(2) = ( r2(1))
      n3(1) = (-r3(2))
      n3(2) = ( r3(1))
      !
      ! edges length
      !
      len(ie,1) = dsqrt(r1(1)**2+r1(2)**2)
      len(ie,2) = dsqrt(r2(1)**2+r2(2)**2)
      len(ie,3) = dsqrt(r3(1)**2+r3(2)**2)
      !
      ! inward normal used for propagation (not normalized)
      !
      ien(ie,1) = n1(1)
      ien(ie,2) = n1(2)
      ien(ie,3) = n2(1)
      ien(ie,4) = n2(2)
      ien(ie,5) = n3(1)
      ien(ie,6) = n3(2)
 
    END DO
 

References fix_periodcity(), w3gdatmd::ien, w3gdatmd::len, w3gdatmd::ntri, w3servmd::strace(), w3gdatmd::trigp, w3gdatmd::xgrd, and w3gdatmd::ygrd.

Referenced by readmsh(), and w3initmd::w3init().

readmsh()

subroutine w3triamd::readmsh	(	integer, intent(in)	NDS,
		character(60), intent(in)	FNAME
	)

Reads triangle and unstructured grid information from GMSH files.

Calls the subroutines needed to compute grid connectivity. Look for namelist with name NAME in unit NDS and read if found.

Parameters

[in]	NDS	Data set number used for search.
[in]	FNAME	Name of namelist.

Author

F. Ardhuin

A. Roland

Date

06-Jun-2018

Definition at line 134 of file w3triamd.F90.

    !/ -------------------------------------------------------------------
    !/                  +-----------------------------------+
    !/                  | WAVEWATCH III           NOAA/NCEP |
    !/                  |           F. Ardhuin              |
    !/                  |           A. Roland               |
    !/                  |                        FORTRAN 90 |
    !/                  | Last update :          06-Jun-2018|
    !/                  +-----------------------------------+
    !/
    !/    15-Feb-2008 : Origination.                        ( version 3.13 )
    !/    25-Aug-2011 : Change of method for IOBPD          ( version 4.04 )
    !/    06-Jun-2018 : Add DEBUGINIT/PDLIB/DEBUGSTP/DEBUGSETIOBP
    !/                                                      ( version 6.04 )
    !/
    !
    !  1. Purpose :
    !
    !      Reads triangle and unstructured grid information from GMSH files
    !      Calls the subroutines needed to compute grid connectivity
    !
    !  2. Method :
    !
    !     Look for namelist with name NAME in unit NDS and read if found.
    !
    !  3. Parameters :
    !
    !     Parameter list
    !     ----------------------------------------------------------------
    !       NDS     Int.   I   Data set number used for search.
    !       NAME    C*4    I   Name of namelist.
    !       STATUS  C*20   O   Status at end of routine,
    !                            '(default values)  ' if no namelist found.
    !                            '(user def. values)' if namelist read.
    !     ----------------------------------------------------------------
    !
    !  4. Subroutines used :
    !
    !      Name               Type  Module   Description
    !     ------------------------------------------------------------------------------------
    !      NEXTLN             Subr.
    !      COUNT              Subr. Internal Count connection.
    !      SPATIAL_GRID       Subr.   Id.    Calculate surfaces.
    !      NVECTRI            Subr.   Id.    Define cell normals and angles and edge length
    !      COORDMAX           Subr.   Id.    Calculate  useful grid elements
    !      AREA_SI            Subr.   Id.    Define Connections
    !     ----------------------------------------------------------------
    !
    !
    !
    !  5. Called by :
    !      Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !      W3GRID    Prog.          Model configuration program
    !     ----------------------------------------------------------------
    !
    !  6. Error messages :
    !
    !  7. Remarks :
    !     The only point index which is needed is IX and NX stands for the total number of grid point.
    !     IY and NY are not needed anymore, they are set to 1 in the unstructured case
    !     Some noticeable arrays are:
    !                     TRIGP  : give the vertices of each triangle
    !     GMSH file gives too much information that is not necessarily required so data processing is needed (data sort and nesting).
    !  8. Structure :
    !
    !  9. Switches :
    !
    ! 10. Source code :
    !
    !/ ------------------------------------------------------------------- /
    USE w3odatmd, ONLY: ndse, ndst, ndso
    USE w3gdatmd, ONLY: zb, xgrd, ygrd, ntri, nx, countot, trigp, nnz, w3dimug
    USE w3servmd, ONLY: itrace, nextln, extcde
    USE constants, only: lpdlib
    USE w3odatmd, ONLY: iaproc
    !
    IMPLICIT NONE
    !/
    !/ Parameter list
    !/
    INTEGER, INTENT(IN)                :: NDS
    CHARACTER(60), INTENT(IN)          :: FNAME
    !/
    !/ local parameters
    !/
    INTEGER                            :: i,j,k, NODES, NELTS, ID, KID
    INTEGER                            :: ID1, ID2, KID1, ITMP(3)
    INTEGER                            :: I1, I2, I3
    INTEGER(KIND=4)                    :: Ind,eltype,ntag, INode
    CHARACTER                          :: COMSTR*1, SPACE*1 = ' ', cels*64
    REAL, ALLOCATABLE                  :: TAGS(:)
    CHARACTER(LEN=64), ALLOCATABLE     :: ELS(:)
    CHARACTER(LEN=120)                 :: LINE
    CHARACTER(LEN=50)                  :: CHTMP
    CHARACTER(LEN=10)                  :: A, B, C
    INTEGER,ALLOCATABLE                :: NELS(:), TRIGPTMP1(:,:), TRIGPTMP2(:,:)
    INTEGER(KIND=4),ALLOCATABLE        :: IFOUND(:), VERTEX(:), BOUNDTMP(:)
    DOUBLE PRECISION, ALLOCATABLE      :: XYBTMP1(:,:),XYBTMP2(:,:)
    REAL                               :: z
 
    OPEN(nds,file = fname,status='old')
    READ (nds,'(A)') comstr
    IF (comstr.EQ.' ') comstr = '$'
    CALL nextln(comstr, nds, ndse)
    READ(nds,*) i,j,k
    CALL nextln(comstr, nds, ndse)
    lpdlib = .false.
#ifdef W3_PDLIB
    lpdlib = .true.
#endif
    !
    ! read number of nodes and nodes from Gmsh files
    !
    READ(nds,*) nodes
    ALLOCATE(xybtmp1(3,nodes))
    DO i= 1, nodes
      READ(nds,*) j, xybtmp1(1,i), xybtmp1(2,i), xybtmp1(3,i)
    END DO
    !
    ! read number of elements and elements from Gmsh files
    !
    ALLOCATE(boundtmp(nodes))
    n_outside_boundary = 0
    CALL nextln(comstr, nds, ndse)
    READ(nds,*) nelts
    ALLOCATE(trigptmp1(3,nelts))
    j = 0
    DO i= 1, nelts
      READ(nds,'(A100)') line
      READ(line,*) ind,eltype,ntag
      ALLOCATE(tags(ntag))
      SELECT CASE (eltype)
        !
        ! eltype = 15 : boundary points  (this is used to make the difference
        !                                between the outside polygon and islands)
        !
      CASE(15)
        READ(line,*) ind,eltype,ntag,tags,inode
        n_outside_boundary = n_outside_boundary +1
        boundtmp(n_outside_boundary)=inode
        !
        ! eltype = 2 : triangles
        !
      CASE (2)
        j = j + 1
        READ(line,*)  ind,eltype,ntag,tags,itmp
        trigptmp1(1:3,j) = itmp
      END SELECT
 
      DEALLOCATE(tags)
    END DO
    !
    ! organizes the grid data structure
    !
    ALLOCATE(outside_boundary(n_outside_boundary))
    outside_boundary(:) = boundtmp(1:n_outside_boundary)
    ntri = j
 
    ALLOCATE(ifound(nodes))
 
    ifound = 0
    !
    ! Verifies that the nodes are used in at least one triangle
    !
    DO k = 1, ntri
      i1 = trigptmp1(1,k)
      i2 = trigptmp1(2,k)
      i3 = trigptmp1(3,k)
 
      ifound(i1)= ifound(i1) + 1
      ifound(i2)= ifound(i2) + 1
      ifound(i3)= ifound(i3) + 1
    END DO
 
    j = 0
 
    ALLOCATE(trigptmp2(3,ntri),vertex(nodes),xybtmp2(3,nodes))
    vertex(:)=0
    xybtmp2 = 0
 
    DO i = 1, nodes
      IF( ifound(i) .GT. 0) THEN
        j = j+1
        xybtmp2(:,j) = xybtmp1(:,i)
        vertex(i) = j
      END IF
    END DO
    !
    ! Number of nodes after clean up
    !
    nx = j
    !
    DO i = 1, ntri
      i1 = trigptmp1(1,i)
      i2 = trigptmp1(2,i)
      i3 = trigptmp1(3,i)
      trigptmp2(1,i) = vertex(i1)
      trigptmp2(2,i) = vertex(i2)
      trigptmp2(3,i) = vertex(i3)
    END DO
    !
    DEALLOCATE(xybtmp1,ifound,trigptmp1)
    DEALLOCATE(vertex)
    !
    !count points connections to allocate array in W3DIMUG
    !
    CALL count(trigptmp2)
    CALL w3dimug ( 1, ntri, nx, countot, nnz, ndse, ndst )
    !
    ! fills arrays
    !
    DO i = 1, nx
      xgrd(1,i) = xybtmp2(1,i)
      ygrd(1,i) = xybtmp2(2,i)
      zb(i)     = xybtmp2(3,i)
    END DO
    !
    DO i=1, ntri
      itmp = trigptmp2(:,i)
      trigp(:,i) = itmp
    END DO
    !
    DEALLOCATE(trigptmp2,xybtmp2)
    !
    ! call the various routines which define the point spotting strategy
    !
    CALL spatial_grid
    CALL nvectri
    CALL coordmax
 
#ifdef W3_PDLIB
    IF(.false.) THEN
#endif
      CALL area_si(1)
#ifdef W3_PDLIB
    ENDIF
#endif
    !
    CLOSE(nds)

References area_si(), coordmax(), count(), w3gdatmd::countot, w3servmd::extcde(), file(), w3odatmd::iaproc, w3servmd::itrace(), constants::lpdlib, n_outside_boundary, w3odatmd::ndse, w3odatmd::ndso, w3odatmd::ndst, w3servmd::nextln(), w3gdatmd::nnz, w3gdatmd::ntri, nvectri(), w3gdatmd::nx, outside_boundary, spatial_grid(), w3gdatmd::trigp, w3gdatmd::w3dimug(), w3gdatmd::xgrd, w3gdatmd::ygrd, and w3gdatmd::zb.

readmsh_iobp()

subroutine w3triamd::readmsh_iobp	(	integer, intent(in)	NDS,
		character(60), intent(in)	FNAME
	)

Reads triangle and unstructured grid information from GMSH files.

Calls the subroutines needed to compute grid connectivity. Look for namelist with name NAME in unit NDS and read if found.

Parameters

[in]	NDS	Data set number used for search.
[in]	FNAME	Name of namelist.

Author

Date

Definition at line 389 of file w3triamd.F90.

    !/ -------------------------------------------------------------------
    !/                  +-----------------------------------+
    !/                  | WAVEWATCH III           NOAA/NCEP |
    !/                  |           A. Roland               |
    !/                  |           F. Ardhuin              |
    !/                  |                        FORTRAN 90 |
    !/                  | Last update :          06-Jun-2018|
    !/                  +-----------------------------------+
    !/
    !/    15-Feb-2008 : Origination.                        ( version 3.13 )
    !/    25-Aug-2011 : Change of method for IOBPD          ( version 4.04 )
    !/    06-Jun-2018 : Add DEBUGINIT/PDLIB/DEBUGSTP/DEBUGSETIOBP
    !/                                                      ( version 6.04 )
    !/
    !
    !  1. Purpose :
    !
    !      Reads triangle and unstructured grid information from GMSH files
    !      Calls the subroutines needed to compute grid connectivity
    !
    !  2. Method :
    !
    !     Look for namelist with name NAME in unit NDS and read if found.
    !
    !  3. Parameters :
    !
    !     Parameter list
    !     ----------------------------------------------------------------
    !       NDS     Int.   I   Data set number used for search.
    !       NAME    C*4    I   Name of namelist.
    !       STATUS  C*20   O   Status at end of routine,
    !                            '(default values)  ' if no namelist found.
    !                            '(user def. values)' if namelist read.
    !     ----------------------------------------------------------------
    !
    !  4. Subroutines used :
    !
    !      Name               Type  Module   Description
    !     ------------------------------------------------------------------------------------
    !      NEXTLN             Subr.
    !      COUNT              Subr. Internal Count connection.
    !      SPATIAL_GRID       Subr.   Id.    Calculate surfaces.
    !      NVECTRI            Subr.   Id.    Define cell normals and angles and edge length
    !      COORDMAX           Subr.   Id.    Calculate  useful grid elements
    !      AREA_SI            Subr.   Id.    Define Connections
    !     ----------------------------------------------------------------
    !
    !
    !
    !  5. Called by :
    !      Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !      W3GRID    Prog.          Model configuration program
    !     ----------------------------------------------------------------
    !
    !  6. Error messages :
    !
    !  7. Remarks :
    !     The only point index which is needed is IX and NX stands for the total number of grid point.
    !     IY and NY are not needed anymore, they are set to 1 in the unstructured case
    !     Some noticeable arrays are:
    !                     TRIGP  : give the vertices of each triangle
    !     GMSH file gives too much information that is not necessarily required so data processing is needed (data sort and nesting).
    !  8. Structure :
    !
    !  9. Switches :
    !
    ! 10. Source code :
    !
    !/ ------------------------------------------------------------------- /
    USE w3odatmd, ONLY: ndse, ndst, ndso
    USE w3gdatmd
    USE w3servmd, ONLY: itrace, nextln, extcde
    USE constants, only: lpdlib
    USE w3odatmd, ONLY: iaproc
    !
    IMPLICIT NONE
    !/
    !/ Parameter list
    !/
    INTEGER, INTENT(IN)                :: NDS
    CHARACTER(60), INTENT(IN)          :: FNAME
    !/
    !/ local parameters
    !/
    INTEGER                            :: i,j,k, NODES
    LOGICAL                            :: lfile_exists
    CHARACTER                          :: COMSTR*1, SPACE*1 = ' ', cels*64
    DOUBLE PRECISION, ALLOCATABLE      :: XYBTMP1(:,:)
 
    INQUIRE(file=fname, exist=lfile_exists)
    IF (.NOT. lfile_exists) RETURN
    OPEN(nds,file = fname,status='old')
    READ (nds,'(A)') comstr
    IF (comstr.EQ.' ') comstr = '$'
    CALL nextln(comstr, nds, ndse)
    READ(nds,*) i,j,k
    CALL nextln(comstr, nds, ndse)
    !
    ! read number of nodes and nodes from Gmsh files
    !
    READ(nds,*) nodes
    ALLOCATE(xybtmp1(3,nodes))
    DO i= 1, nodes
      READ(nds,*) j, xybtmp1(1,i), xybtmp1(2,i), xybtmp1(3,i)
      IF (int(xybtmp1(3,i)) .EQ. 3) iobp(i) = 3
    END DO
    !
    CLOSE(nds)

References w3servmd::extcde(), file(), w3odatmd::iaproc, w3gdatmd::iobp, w3servmd::itrace(), constants::lpdlib, w3odatmd::ndse, w3odatmd::ndso, w3odatmd::ndst, and w3servmd::nextln().

Referenced by set_ug_iobp().

readmshobc()

subroutine w3triamd::readmshobc	(	integer, intent(in)	NDS,
		character(60), intent(in)	FNAME,
		integer, dimension(ny,nx), intent(inout)	TMPSTA,
		logical, intent(out)	UGOBCOK
	)

Reads open boundary information for UNST grids following GMESH type format.

Parameters

[in]	NDS	Data set number used for search.
[in]	FNAME	File name.
[in,out]	TMPSTA	Status map to be updated (for OBC,TMPSTA = 2).
[out]	UGOBCOK	Flag for proper reading of OBC file.

Author

F. Ardhuin

Date

14-Mar-2018

Definition at line 670 of file w3triamd.F90.

    !/ -------------------------------------------------------------------
    !/                  +-----------------------------------+
    !/                  | WAVEWATCH III           NOAA/NCEP |
    !/                  |           F. Ardhuin              |
    !/                  |                        FORTRAN 90 |
    !/                  | Last update :          14-Mar-2018|
    !/                  +-----------------------------------+
    !/
    !/    14-Mar-2018 : Origination.                        ( version 6.02 )
    !/
    !
    !  1. Purpose :
    !
    !      Reads open boundary information for UNST grids following GMESH type format
    !
    !  2. Method :
    !
    !     Reads an ASCII file
    !
    !  3. Parameters :
    !
    !     Parameter list
    !     ----------------------------------------------------------------
    !       NDS     Int.   I   Data set number used for search.
    !       FNAME   Char*60 I  File name
    !       TMPSTA  Char*60 I/O status map to be updated (for OBC, TMPSTA = 2)
    !       UGOBCOK Logical O   flag for proper reading of OBC file
    !     ----------------------------------------------------------------
    !
    !  4. Subroutines used : NONE
    !
    !  5. Called by :
    !      Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !      W3GRID    Prog.          Model configuration program
    !     ----------------------------------------------------------------
    !
    !  6. Error messages :
    !
    !  7. Remarks :
    !
    !  8. Structure :
    !
    !  9. Switches :
    !
    ! 10. Source code :
    !
    !/ ------------------------------------------------------------------- /
    USE w3gdatmd, ONLY: nx, ny, ccon , countcon
    USE w3odatmd, ONLY: ndse, ndst, ndso
    USE w3servmd, ONLY: itrace, nextln, extcde
#ifdef W3_S
    USE w3servmd, ONLY: strace
#endif
 
    IMPLICIT NONE
    !/ ------------------------------------------------------------------- /
    !/ Parameter list
    !/
    INTEGER, INTENT(IN)                :: NDS
    CHARACTER(60), INTENT(IN) :: FNAME
    INTEGER, INTENT(INOUT)    :: TMPSTA(NY,NX)
    LOGICAL, INTENT(OUT)      :: UGOBCOK
    !/
    !/ local parameters
    !/
    INTEGER                            :: I, IERR
    INTEGER(KIND=4)                    :: Ind,ntag, INode
    CHARACTER                          :: COMSTR*1, SPACE*1 = ' ', cels*64
    REAL, ALLOCATABLE                  :: TAGS(:)
    CHARACTER(LEN=120)                 :: LINE
 
    ugobcok=.false.
 
    OPEN(nds,file = fname,status='old')
    READ (nds,'(A)') comstr
    IF (comstr.EQ.' ') comstr = '$'
    CALL nextln(comstr, nds, ndse)
    ierr = 0
    DO WHILE (ierr.EQ.0)
      READ (nds,'(A100)',END=2001,ERR=2002,IOSTAT=IERR) line
      READ(line,*,iostat=ierr) ind,ntag
      IF (ierr.EQ.0) THEN
        ALLOCATE(tags(ntag))
        READ(line,*,iostat=ierr) ind,ntag,tags,inode
        IF (ierr.EQ.0) THEN
          tmpsta(1,inode)=2
          DEALLOCATE(tags)
        ELSE
          GOTO 2001
        END IF
      END IF
    END DO
    CLOSE(nds)
    ugobcok=.true.
    RETURN
    !
2001 CONTINUE
    WRITE (ndse,1001)
    CALL extcde ( 61 )
    !
2002 CONTINUE
    WRITE (ndse,1002) ierr
    CALL extcde ( 62 )
1001 FORMAT (/' *** WAVEWATCH III ERROR IN READMSHOBC : '/          &
         '     PREMATURE END OF FILE IN READING ',a/)
1002 FORMAT (/' *** WAVEWATCH III ERROR IN READMSHOBC : '/          &
         '     ERROR IN READING ',a,'  IOSTAT =',i8/)
 

References w3gdatmd::ccon, w3gdatmd::countcon, w3servmd::extcde(), file(), w3servmd::itrace(), w3odatmd::ndse, w3odatmd::ndso, w3odatmd::ndst, w3servmd::nextln(), w3gdatmd::nx, w3gdatmd::ny, and w3servmd::strace().

set_iobp()

subroutine w3triamd::set_iobp	(	integer, dimension(nx), intent(in)	MASK,
		integer*2, dimension(nx), intent(out)	STATUS
	)

Setup boundary pointer.

Parameters

[in,out]	MASK
[in,out]	STATUS

Author

Aron Roland

Mathiew Dutour-Sikiric

Date

01-Jun-2018

Definition at line 2319 of file w3triamd.F90.

    !/
    !/                  +-----------------------------------+
    !/                  | WAVEWATCH III           NOAA/NCEP |
    !/                  |                                   |
    !/                  | Aron Roland (BGS IT&E GmbH)       |
    !/                  | Mathieu Dutour-Sikiric (IRB)      |
    !/                  |                                   |
    !/                  |                        FORTRAN 90 |
    !/                  | Last update :        01-June-2018 |
    !/                  +-----------------------------------+
    !/
    !/    01-June-2018 : Origination.                        ( version 6.04 )
    !/
    !  1. Purpose : setup boundary pointer
    !  2. Method :
    !  3. Parameters :
    !
    !     Parameter list
    !     ----------------------------------------------------------------
    !     ----------------------------------------------------------------
    !
    !  4. Subroutines used :
    !
    !      Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !      STRACE    Subr. W3SERVMD Subroutine tracing.
    !     ----------------------------------------------------------------
    !
    !  5. Called by :
    !
    !      Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !     ----------------------------------------------------------------
    !
    !  6. Error messages :
    !  7. Remarks
    !  8. Structure :
    !  9. Switches :
    !
    !     !/S  Enable subroutine tracing.
    !
    ! 10. Source code :
    !
    !/ ------------------------------------------------------------------- /
#ifdef W3_S
    USE w3servmd, ONLY: strace
#endif
    !
    !/
    !
    USE constants
    !
    !
    USE w3gdatmd, ONLY: nx, ntri, trigp
    USE w3odatmd, ONLY: iaproc
 
 
    IMPLICIT NONE
 
    !/
    !/ ------------------------------------------------------------------- /
    !/ Parameter list
    !/
    INTEGER, INTENT(IN)   :: MASK(NX)
    INTEGER*2, INTENT(OUT)  :: STATUS(NX)
    !
    INTEGER :: COLLECTED(NX), NEXTVERT(NX), PREVVERT(NX)
    INTEGER          :: ISFINISHED !, INEXT, IPREV
    INTEGER :: INEXT(3), IPREV(3)
    INTEGER          :: ZNEXT, IP, I, IE, IPNEXT, IPPREV, COUNT
    integer nb0, nb1, nbM1
    status = -1
    inext=(/ 2, 3, 1 /) !IPREV=1+MOD(I+1,3)
    iprev=(/ 3, 1, 2 /) !INEXT=1+MOD(I,3)
    DO ie=1,ntri
      ! If one of the points of the triangle is masked out (land) then do as if triangle does not exist...
      !        IF ((MASK(TRIGP(1,IE)).GT.0).AND.(MASK(TRIGP(2,IE)).GT.0).AND.(MASK(TRIGP(3,IE)).GT.0)) THEN
      DO i=1,3
        ip=trigp(i,ie)
        CALL triang_indexes(i, ipnext, ipprev)
        !IPNEXT=TRIGP(INEXT(I),IE)
        !IPPREV=TRIGP(IPREV(I),IE)
        IF (status(ip).EQ.-1) THEN
          status(ip)=1
          prevvert(ip)=ipprev
          nextvert(ip)=ipnext
        END IF
      END DO
      !        ENDIF
    END DO
    status(:)=-1
    !
    count = 0
    DO
      count = count + 1
      collected(:)=0
      DO ie=1,ntri
        !        IF ((MASK(TRIGP(1,IE)).GT.0).AND.(MASK(TRIGP(2,IE)).GT.0).AND.(MASK(TRIGP(3,IE)).GT.0)) THEN
        DO i=1,3
          ip=trigp(i,ie)
          CALL triang_indexes(i, ipnext, ipprev)
          !IPNEXT=TRIGP(INEXT(I),IE)
          !IPPREV=TRIGP(IPREV(I),IE)
          IF (status(ip).EQ.-1) THEN
            znext=nextvert(ip)
            IF (znext.EQ.ipprev) THEN
              collected(ip)=1
              nextvert(ip)=ipnext
              IF (nextvert(ip).EQ.prevvert(ip)) THEN
                status(ip)=1
              END IF
            END IF
          END IF
        END DO
        !            END IF ! end of test on MASK
      END DO
      !
      ! Checks that all nodes have been treated ...
      !
      isfinished=1
      DO ip=1,nx
        IF (mask(ip).LE.0) THEN
          status(ip)=0
        ELSE
          IF ((collected(ip).EQ.0).AND.(status(ip).EQ.-1)) THEN
            status(ip)=0
          END IF
          IF (status(ip).eq.-1) THEN
            isfinished=0
          END IF
        ENDIF
      END DO
      IF (isfinished.EQ.1) THEN
        EXIT
      END IF
    END DO
 
    status = 1
    CALL get_boundary(nx, ntri, trigp, status, prevvert, nextvert)
 
    !#ifdef MPI_PARALL_GRID
    !      CALL exchange_p2di(STATUS)
    !#endif

References get_boundary(), w3odatmd::iaproc, w3gdatmd::ntri, w3gdatmd::nx, w3servmd::strace(), triang_indexes(), and w3gdatmd::trigp.

Referenced by set_ug_iobp(), and ug_getopenboundary().

set_ug_iobp()

subroutine w3triamd::set_ug_iobp

Redefines the values of the boundary points and angle pointers based on the MAPSTA array.

Adapted boundary detection from A. Roland and M. Dutour (WWM code).

Author

Fabrice Ardhuin

Aron Roland

Date

17-Apr-2016

Definition at line 2811 of file w3triamd.F90.

    !/
    !/                  +-----------------------------------+
    !/                  | WAVEWATCH III           NOAA/NCEP |
    !/                  |        Fabrice Ardhuin            |
    !/                  |        Aron Roland                |
    !/                  |                        FORTRAN 90 |
    !/                  | Last update :         17-Apr-2016 |
    !/                  +-----------------------------------+
    !/
    !/    23-Aug-2011 : Origination.                        ( version 4.04 )
    !/    17-Apr-2016 : Uses optimized boundary detection   ( version 5.10 )
    !/
    !  1. Purpose :
    !
    !     Redefines the values of the boundary points and angle pointers
    !     based on the MAPSTA array
    !
    !  2. Method :
    !
    !     Adapted boundary detection from A. Roland and M. Dutour (WWM code)
    !
    !  3. Parameters :
    !
    !     Parameter list
    !     ----------------------------------------------------------------
    !     ----------------------------------------------------------------
    !
    !     Local variables.
    !     ----------------------------------------------------------------
    !     ----------------------------------------------------------------
    !
    !  4. Subroutines used :
    !
 
    !  5. Called by :
    !
    !      Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !      WW3_GRID  Prog. WW3_GRID Grid preprocessor
    !      W3ULEV    Subr. W3UPDTMD Water level update
    !     ----------------------------------------------------------------
    !
    !  6. Error messages :
    !
    !       None.
    !
    !  7. Remarks :
    !
    !  8. Structure :
    !
    !
    !  9. Switches :
    !
    !       !/S     Enable subroutine tracing.
    !
    !
    ! 10. Source code :
    !/ ------------------------------------------------------------------- /
    !/
    !
    USE constants
    !
    !
    USE w3gdatmd, ONLY: nx, ny, nsea, mapfs,                        &
         nk, nth, dth, xfr, mapsta, countri,         &
         ecos, esin, ien, ntri, trigp,               &
         iobp,iobpd, iobpa,                          &
#ifdef w3_ref1
         refpars, reflc, refld,                      &
#endif
         angle0, angle
 
    USE w3odatmd, ONLY: tbpi0, tbpin, flbpi
    USE w3adatmd, ONLY: cg, cx, cy, atrnx, atrny, itime, cflxymax
    USE w3idatmd, ONLY: flcur
    USE w3odatmd, only : iaproc
#ifdef W3_S
    USE w3servmd, ONLY: strace
#endif
 
    IMPLICIT NONE
    !/ ------------------------------------------------------------------- /
    !/ Parameter list
    !/
    !/
    !/ ------------------------------------------------------------------- /
    !/ Local parameters
    !/
    INTEGER                 :: ITH, IX, I, J, IP, IE, NDIRSUM
    REAL (KIND = 8)         :: cossum, sinsum
    REAL (KIND = 8)         :: dirmin, dirmax, shift, tempo, dircoast
    REAL (KIND = 8)         :: x1, x2, y1, y2, dxp1, dxp2, dxp3
    REAL (KIND = 8)         :: dyp1, dyp2, dyp3, edet1, edet2, evx, evy
    REAL(KIND=8), parameter :: thr    = tiny(1.)
    INTEGER                 :: I1, I2, I3
    INTEGER                 :: ITMP(NX), NEXTVERT(NX), PREVVERT(NX)
    CHARACTER(60) :: FNAME
#ifdef W3_S
    INTEGER, SAVE           :: IENT = 0
#endif
    !/ ------------------------------------------------------------------- /
    !
    ! 1.  Preparations --------------------------------------------------- *
    ! 1.a Set constants
    !
#ifdef W3_S
    CALL strace (ient, 'SETUGIOBP')
#endif
    !
    !--- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
    ! 2.  Searches for boundary points
    !
    itmp = mapsta(1,:)
    CALL set_iobp(itmp, iobp)
    fname = 'meshbnd.msh'
    CALL readmsh_iobp(23456,fname)
    !
    !--- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
    ! 3. Defines directions pointing into land or sea
    !
    iobpd(:,:) = 0
    iobpa(:)   = 0
    !
    DO ip=1,nx
      IF (mapsta(1,ip).EQ.2) THEN
        iobpa(ip) = 1
        iobp(ip)  = 2
      ENDIF
    END DO
 
    DO ie = 1,ntri
      i1   =   trigp(1,ie)
      i2   =   trigp(2,ie)
      i3   =   trigp(3,ie)
      dxp1 =   ien(ie,6)
      dyp1 = - ien(ie,5)
      dxp2 =   ien(ie,2)
      dyp2 = - ien(ie,1)
      dxp3 =   ien(ie,4)
      dyp3 = - ien(ie,3)
      DO ith=1,nth
        evx=ecos(ith)
        evy=esin(ith)
        DO i=1,3
          IF (i.eq.1) THEN
            x1=   dxp1
            y1=   dyp1
            x2= - dxp3
            y2= - dyp3
            ip=   i1
          END IF
          IF (i.eq.2) THEN
            x1 =   dxp2
            y1 =   dyp2
            x2 = - dxp1
            y2 = - dyp1
            ip =   i2
          END IF
          IF (i.eq.3) THEN
            x1 =   dxp3
            y1 =   dyp3
            x2 = - dxp2
            y2 = - dyp2
            ip =   i3
          END IF
          IF (iobp(ip) .eq. 0) THEN ! physical boundary
            edet1 = thr-x1*evy+y1*evx
            edet2 = thr+x2*evy-y2*evx
            IF ((edet1.gt.0.).and.(edet2.gt.0.)) THEN
              ! this is the case of waves going towards the boundary ...
              iobpd(ith,ip)=1
            ENDIF
          ELSE ! water ...
            iobpd(ith,ip)=1
          END IF
        END DO
      END DO
    END DO
    DO ip = 1, nx
      IF ( iobpa(ip) .eq. 1 .OR. iobp(ip) .eq. 3 .OR. iobp(ip) .eq. 4) iobpd(:,ip) = 1
    END DO
    !2do: recode for mpi
    !        IF (LBCWA .OR. LBCSP) THEN
    !          IF (.NOT. ANY(IOBP .EQ. 2)) THEN
    !            CALL WWM_ABORT('YOU IMPOSED BOUNDARY CONDITIONS BUT IN THE BOUNDARY FILE ARE NO NODES WITH FLAG = 2')
    !          ENDIF
    !        ENDIF
    !#ifdef MPI_PARALL_GRID
    !      CALL exchange_p2di(IOBWB)
    !      DO ID = 1, MDC
    !        iwild = IOBPD(ID,:)
    !        CALL exchange_p2di(iwild)
    !        IOBPD(ID,:) = iwild
    !      ENDDO
    !#endif
    !--- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
    ! 3. Updates the reflection direction and sharp / flat shoreline angle
 
#ifdef W3_REF1
    !
    ! Finds the shoreline direction from IOBPD
    !
    reflc(1,:)= 0.
    refld(:,:)= 1
    DO ip=1,nx
      IF (iobp(ip).EQ.0.AND.mapsta(1,ip).EQ.1) THEN
        cossum=0.
        sinsum=0.
        ndirsum=0.
        DO ith=1,nth
          cossum=cossum+iobpd(ith,ip)*ecos(ith)
          sinsum=sinsum+iobpd(ith,ip)*esin(ith)
          ndirsum=ndirsum+iobpd(ith,ip)
        END DO
        dircoast=atan2(sinsum, cossum)
        refld(1,mapfs(1,ip)) = 1+mod(nth+nint(dircoast/dth),nth)
        refld(2,mapfs(1,ip)) = 4-max(2,nint(4.*real(ndirsum)/real(nth)))
        reflc(1,mapfs(1,ip))= refpars(1)
      END IF
    END DO
#endif
    !
    ! Recomputes the angles used in the gradients estimation
    !
    !
    RETURN

References w3adatmd::atrnx, w3adatmd::atrny, w3adatmd::cflxymax, w3adatmd::cg, w3gdatmd::countri, w3adatmd::cx, w3adatmd::cy, w3gdatmd::dth, w3gdatmd::ecos, w3gdatmd::esin, w3odatmd::flbpi, w3idatmd::flcur, w3odatmd::iaproc, w3gdatmd::ien, w3gdatmd::iobp, w3gdatmd::iobpa, w3gdatmd::iobpd, w3adatmd::itime, w3gdatmd::mapfs, w3gdatmd::mapsta, w3gdatmd::nk, w3gdatmd::nsea, w3gdatmd::nth, w3gdatmd::ntri, w3gdatmd::nx, w3gdatmd::ny, readmsh_iobp(), set_iobp(), w3servmd::strace(), w3odatmd::tbpi0, w3odatmd::tbpin, w3gdatmd::trigp, and w3gdatmd::xfr.

Referenced by w3iorsmd::w3iors().

spatial_grid()

subroutine w3triamd::spatial_grid

Calculates triangle areas and reorders the triangles to have them oriented counterclockwise.

The triangle surface calculation is based on cross product.

Author

A. Roland

F. Ardhuin

Date

31-Aug-2011

Definition at line 891 of file w3triamd.F90.

    !/ -------------------------------------------------------------------
    !/                  +-----------------------------------+
    !/                  | WAVEWATCH III           NOAA/NCEP |
    !/                  |      A. Roland  (BGS IT&E GbmH)   |
    !/                  |      F. Ardhuin (IFREMER)         |
    !/                  |                        FORTRAN 90 |
    !/                  | Last update :          31-Aug-2011|
    !/                  +-----------------------------------+
    !/
    !/    15-May-2007 : Origination: adjustment from the WWM code       ( version 3.13 )
    !/    31-Aug-2011 : Simplfies the cross products                    ( version 4.05 )
    !/
    !
    !  1. Purpose :
    !
    !      Calculates triangle areas and reorders the triangles to have them
    !      oriented counterclockwise
    !
    !  2. Method :
    !
    !     The triangle surface calculation is based on cross product.
    !
    !  3. Parameters :
    !
    !  4. Subroutines used :
    !
    !  5. Called by :
    !
    !       Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !      READTRI    Subr. Internal  Unstructured mesh definition.
    !     ----------------------------------------------------------------
    !
    !  6. Error messages :
    !
    !  7. Remarks :
    !
    !     This part of code is adapted from the WWM wave model develop at the Darmstadt University
    !     (Aaron Roland)
    !
    !  8. Structure :
    !
    !  9. Switches :
    !
    ! 10. Source code :
    !
    !/ ------------------------------------------------------------------- /
    USE w3gdatmd
#ifdef W3_S
    USE w3servmd, ONLY: strace
#endif
    USE w3odatmd, ONLY: ndse
 
    IMPLICIT NONE
    !
    !local parameters
    !
    REAL              :: TL1, TL2, TL3, TMPTRIGP
    INTEGER           :: I1, I2, I3
    INTEGER           :: K
    real*8            :: pt(3,2)
 
#ifdef W3_S
    INTEGER                      ::  IENT = 0
#endif
    !/ ------------------------------------------------------------------- /
#ifdef W3_S
    CALL strace (ient, 'SPATIAL_GRID')
#endif
 
    DO k = 1, ntri
 
      i1 = trigp(1,k)
      i2 = trigp(2,k)
      i3 = trigp(3,k)
 
!AR: todo call this only for global grid 
      CALL fix_periodcity(i1,i2,i3,xgrd,ygrd,pt)
      !
      ! cross product of edge-vector  (orientated anticlockwise)
      !
      tria(k) = real( (pt(2,2)-pt(1,2)) & 
           *(pt(1,1)-pt(3,1))      &    
           +(pt(3,2)-pt(1,2))      &    
           *(pt(2,1)-pt(1,1))      )*0.5
      !
      ! test on negative triangle area, which means that the orientiation is not as assumed to be anticw.
      ! therefore we swap the nodes !!!
      !
      IF (tria(k) .lt. tiny(1.)) THEN
        tmptrigp = trigp(2,k)
        trigp(2,k) = trigp(3,k)
        trigp(3,k) = tmptrigp
        i2 = trigp(2,k)
        i3 = trigp(3,k)
        tria(k) = -1.d0*tria(k)
      END IF
    END DO

References fix_periodcity(), w3odatmd::ndse, w3gdatmd::ntri, w3servmd::strace(), w3gdatmd::tria, w3gdatmd::trigp, w3gdatmd::xgrd, and w3gdatmd::ygrd.

Referenced by readmsh(), and w3initmd::w3init().

triang_indexes()

subroutine w3triamd::triang_indexes	(	integer, intent(in)	I,
		integer, intent(out)	INEXT,
		integer, intent(out)	IPREV
	)

Set indices of the triangle.

Parameters

[in]	I
[out]	INEXT
[out]	IPREV

Author

Aron Roland

Mathieu Dutour-Sikiric

Date

01-Jun-2018

Definition at line 2720 of file w3triamd.F90.

    !/
    !/                  +-----------------------------------+
    !/                  | WAVEWATCH III           NOAA/NCEP |
    !/                  |                                   |
    !/                  | Aron Roland (BGS IT&E GmbH)       |
    !/                  | Mathieu Dutour-Sikiric (IRB)      |
    !/                  |                                   |
    !/                  |                        FORTRAN 90 |
    !/                  | Last update :        01-June-2018 |
    !/                  +-----------------------------------+
    !/
    !/    01-June-2018 : Origination.                        ( version 6.04 )
    !/
    !  1. Purpose : set indices of the triangle
    !  2. Method :
    !  3. Parameters :
    !
    !     Parameter list
    !     ----------------------------------------------------------------
    !     ----------------------------------------------------------------
    !
    !  4. Subroutines used :
    !
    !      Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !      STRACE    Subr. W3SERVMD Subroutine tracing.
    !     ----------------------------------------------------------------
    !
    !  5. Called by :
    !
    !      Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !     ----------------------------------------------------------------
    !
    !  6. Error messages :
    !  7. Remarks
    !  8. Structure :
    !  9. Switches :
    !
    !     !/S  Enable subroutine tracing.
    !
    ! 10. Source code :
    !
    !/ ------------------------------------------------------------------- /
#ifdef W3_S
    USE w3servmd, ONLY: strace
#endif
    IMPLICIT NONE
    !
    !/ ------------------------------------------------------------------- /
    !/ Parameter list
    !/
    !/ ------------------------------------------------------------------- /
    !/ Local PARAMETERs
    !/
#ifdef W3_S
    INTEGER, SAVE           :: IENT = 0
#endif
    !/
    !/ ------------------------------------------------------------------- /
    !/
    INTEGER, INTENT(IN)  :: I
    INTEGER, INTENT(OUT) :: INEXT, IPREV
#ifdef W3_S
    CALL strace (ient, 'TRIANG_INDEXES')
#endif
    IF (i.EQ.1) THEN
      inext=3
    ELSE
      inext=i-1
    END IF
    IF (i.EQ.3) THEN
      iprev=1
    ELSE
      iprev=i+1
    END IF

References w3servmd::strace().

Referenced by get_boundary(), and set_iobp().

ug_getopenboundary()

subroutine w3triamd::ug_getopenboundary	(	integer, dimension(ny,nx), intent(inout)	TMPSTA,
		real, dimension(ny,nx), intent(in)	ZBIN,
		real, intent(in)	ZLIM
	)

Defines open boundary points based on depth.

A boundary node has more node around it than triangles.

Parameters

[in,out]	TMPSTA	Status map to be updated (for OBC, TMPSTA = 2).
[in]	ZBIN
[in]	ZLIM

Author

F. Ardhuin

Date

30-Aug-2012

Definition at line 797 of file w3triamd.F90.

    !/                  +-----------------------------------+
    !/                  | WAVEWATCH III           NOAA/NCEP |
    !/                  |           F. Ardhuin              |
    !/                  |                        FORTRAN 90 |
    !/                  | Last update :          30-Aug-2012|
    !/                  +-----------------------------------+
    !/
    !/    30-Aug-2012 : Adpatation from SHOM-Ifremer program( version 4.07 )
    !/
    !
    !  1. purpose: defines open boundary points based on depth
    !  2. Method : a boundary node has more node around it than triangles
    !
    !
    !  3. Parameters :
    !     TMPSTA: status map to be updated (for OBC, TMPSTA = 2)
    !
    !  4. Subroutines used :
    !
    !  5. Called by :
    !
    !       Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !      w3GRID    Prog.          Model configuration program
    !     ----------------------------------------------------------------
    !
    !  6. Error messages :
    !
    !  7. Remarks :
    !
 
    !
    !  8. Structure :
    !
    !  9. Switches :
    !
    ! 10. Source code :
    USE w3gdatmd, ONLY: nx, ny, ccon, countcon, iobp
 
#ifdef W3_S
    USE w3servmd, ONLY: strace
#endif
 
    IMPLICIT NONE
    !/ ------------------------------------------------------------------- /
    !/ Parameter list
    !/
    INTEGER, INTENT(INOUT)  :: TMPSTA(NY,NX)
#ifdef W3_S
    INTEGER, SAVE           :: IENT = 0
#endif
    REAL   , INTENT(IN)     :: ZBIN(NY,NX)
    REAL   , INTENT(IN)     :: ZLIM
    !/
    !/ ------------------------------------------------------------------- /
    !/ Local parameters
    !/
    INTEGER                 ::  IBC, IX
    INTEGER                 ::  MASK(NX)
    INTEGER*2               ::  STATUS(NX)
    !
    mask(:)=1
    CALL set_iobp (mask, status)
    !
#ifdef W3_S
    CALL strace (ient, 'UG_GETOPENBOUNDARY')
#endif
    DO ibc = 1, n_outside_boundary
      ix = outside_boundary(ibc)
      !write(*,*) 'TEST1', IX, TMPSTA(1,IX), CCON(IX), COUNTCON(IX), ZBIN(1,IX), ZLIM
      ! OUTSIDE_BOUNDARY(IBC) is defined over the full nodes NODES indexes
      ! whereas TMPSTA and ZBIN are defined over the clean up list of nodes NX
      IF ((ix.NE.0).AND.(ix.LE.nx)) THEN
        IF ( (tmpsta(1,ix).EQ.1) .AND. (status(ix).EQ.0) .AND. (zbin(1,ix) .LT. zlim)) tmpsta(1,ix) = 2
      END IF
    END DO
    !

References w3gdatmd::ccon, w3gdatmd::countcon, w3gdatmd::iobp, n_outside_boundary, w3gdatmd::nx, w3gdatmd::ny, outside_boundary, set_iobp(), and w3servmd::strace().

ug_gradients()

subroutine w3triamd::ug_gradients	(	real, dimension(0:nsea), intent(in)	PARAM,
		real, dimension(:,:), intent(out)	DIFFX,
		real, dimension(:,:), intent(out)	DIFFY
	)

Calculate gradients at a point via its connection.

Using linear shape function this is a basis on which all advection schemes in Roland (2008) are checked.

Parameters

[in]	PARAM	Depth or current field (indices 0 to NSEA).
[out]	DIFFX	X gradient (indices 1 to NX).
[out]	DIFFY	Y gradient (indices 1 to NY).

Author

F. Ardhuin

A. Roland

Date

14-Oct-2013

Definition at line 2039 of file w3triamd.F90.

    !/                  +-----------------------------------+
    !/                  | WAVEWATCH III           NOAA/NCEP |
    !/                  |           F. Ardhuin              |
    !/                  |           A. Roland               |
    !/                  |                        FORTRAN 90 |
    !/                  | Last update :          14-Oct-2013|
    !/                  +-----------------------------------+
    !/
    !/    15-Nov-2007 : Origination.                        ( version 3.13 )
    !/    31-Oct-2010 : Merging of 4.03 with 3.14-Ifremer   ( version 4.04 )
    !/    08-Nov-2011 : Correction for zero grad. on contour( version 4.04 )
    !/    14-Oct-2013 : Correction  of latitude factor      ( version 4.12 )
    !/    01-Mai-2018 : Using linear shape function for gradients [ version 6.04)
    !/
    !
    !  1. purpose: calculate gradients at a point via its connection.
    !  2. Method : using linear shape function this is a basis on which
    !              all advection schemes in Roland (2008) are checked.
    !
    !  3. Parameters :
    !     PARAM : depth or current field (indices 0 to NSEA)
    !     DIFFX :  x gradient            (indices 1 to NX)
    !     DIFFY :  y gradient            (indices 1 to NX)
    !
    !  4. Subroutines used :
    !
    !  5. Called by :
    !
    !       Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !      W3WAVE    Subr.          Actual wind wave routine
    !     ----------------------------------------------------------------
    !
    !  6. Error messages :
    !
    !  7. Remarks :
    !
    !      This subroutine is adjusted from WWM code (Aaron Roland)
    !
    !  8. Structure :
    !
    !  9. Switches :
    !
    ! 10. Source code :
    USE constants
    USE w3gdatmd, ONLY : trigp, ntri, nx, nsea, mapfs, clatis, &
         mapsta, angle, flagll,  iobp, ien, tria, nseal, ntri
    USE w3adatmd, ONLY : nsealm
#ifdef W3_PDLIB
    USE yowelementpool
    use yownodepool,    only: pdlib_ien, pdlib_tria, npa
    USE yowexchangemodule, only : pdlib_exchange1dreal
#endif
 
    IMPLICIT NONE
 
 
    REAL, INTENT(IN)     :: PARAM(0:NSEA)
    REAL, INTENT(OUT)  :: DIFFX(:,:), DIFFY(:,:)
 
    ! local parameters
 
    INTEGER              :: VERTICES(3), NI(3), NI_GL(3)
    REAL                 :: TMP1(3), TMP2(3)
    INTEGER              :: I, IX, IE, IE_GL
    REAL                 :: VAR(3), FACT, LATMEAN
    REAL                 :: DIFFXTMP, DIFFYTMP
    REAL                 :: DEDX(3), DEDY(3)
    REAL                 :: DVDXIE, DVDYIE
    REAL                 :: WEI(NX), WEI_LOCAL(NSEAL)
    real*8               :: rtmp(nseal)
 
    diffx = 0.
    diffy = 0.
    !
    IF (flagll) THEN
      fact=1./(dera*radius)
    ELSE
      fact=1.
    END IF
 
#ifdef W3_PDLIB
    IF (.NOT. lpdlib) THEN
#endif
      wei = 0.
      DO ie = 1, ntri
        ni = trigp(:,ie)
        latmean = 1./3. * sum(clatis(mapfs(1,ni)))
        wei(ni) = wei(ni) + 2.*tria(ie)
        dedx(1)     = ien(ie,1)
        dedx(2)     = ien(ie,3)
        dedx(3)     = ien(ie,5)
        dedy(1)     = ien(ie,2)
        dedy(2)     = ien(ie,4)
        dedy(3)     = ien(ie,6)
        var         = param(mapfs(1,ni)) * fact
        dvdxie      = dot_product( var,dedx)
        dvdyie      = dot_product( var,dedy)
        diffx(1,ni) = diffx(1,ni) + dvdxie * latmean
        diffy(1,ni) = diffy(1,ni) + dvdyie
      END DO
      diffx(1,:) = diffx(1,:)/wei
      diffy(1,:) = diffy(1,:)/wei
#ifdef W3_PDLIB
    ELSE
      wei_local = 0.
      DO ie = 1, ne
        ni      = ine(:,ie)
        ie_gl   = ielg(ie)
        ni_gl   = trigp(:,ie_gl)
        latmean = 1./3. * sum(clatis(mapfs(1,ni_gl)))
        wei_local(ni) = wei_local(ni) + 2.*pdlib_tria(ie)
        dedx(1)     = pdlib_ien(1,ie)
        dedx(2)     = pdlib_ien(3,ie)
        dedx(3)     = pdlib_ien(5,ie)
        dedy(1)     = pdlib_ien(2,ie)
        dedy(2)     = pdlib_ien(4,ie)
        dedy(3)     = pdlib_ien(6,ie)
        var         = param(mapfs(1,ni_gl)) * fact
        dvdxie      = dot_product(var,dedx)
        dvdyie      = dot_product(var,dedy)
        diffx(1,ni) = diffx(1,ni) + dvdxie * latmean
        diffy(1,ni) = diffy(1,ni) + dvdyie
      END DO
      diffx(1,:) = diffx(1,:)/wei_local
      diffy(1,:) = diffy(1,:)/wei_local
    ENDIF
    CALL pdlib_exchange1dreal(diffx(1,:))
    CALL pdlib_exchange1dreal(diffy(1,:))
#endif
    !

References w3gdatmd::angle, w3gdatmd::clatis, constants::dera, w3gdatmd::flagll, yowelementpool::ielg, w3gdatmd::ien, yowelementpool::ine, w3gdatmd::iobp, constants::lpdlib, w3gdatmd::mapfs, w3gdatmd::mapsta, yowelementpool::ne, yownodepool::npa, w3gdatmd::nsea, w3gdatmd::nseal, w3adatmd::nsealm, w3gdatmd::ntri, w3gdatmd::nx, yowexchangemodule::pdlib_exchange1dreal(), yownodepool::pdlib_ien, yownodepool::pdlib_tria, constants::radius, w3gdatmd::tria, and w3gdatmd::trigp.

Referenced by w3wavemd::w3wave().

w3nestug()

subroutine w3triamd::w3nestug	(	real, intent(in)	DISTMIN,
		logical, intent(inout)	FLOK
	)

UGTYPE nesting initialization.

Parameters

[in]	DISTMIN
[in,out]	FLOK

Author

Aron Roland

Mathieu Dutour-Sikiric

Date

01-Jun-2018

Definition at line 2184 of file w3triamd.F90.

    !/
    !/                  +-----------------------------------+
    !/                  | WAVEWATCH III           NOAA/NCEP |
    !/                  |                                   |
    !/                  | Aron Roland (BGS IT&E GmbH)       |
    !/                  | Mathieu Dutour-Sikiric (IRB)      |
    !/                  |                                   |
    !/                  |                        FORTRAN 90 |
    !/                  | Last update :        01-June-2018 |
    !/                  +-----------------------------------+
    !/
    !/    01-June-2018 : Origination.                        ( version 6.04 )
    !/
    !  1. Purpose : UGTYPE nesting initialization
    !  2. Method :
    !  3. Parameters :
    !
    !     Parameter list
    !     ----------------------------------------------------------------
    !     ----------------------------------------------------------------
    !
    !  4. Subroutines used :
    !
    !      Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !      STRACE    Subr. W3SERVMD Subroutine tracing.
    !     ----------------------------------------------------------------
    !
    !  5. Called by :
    !
    !      Name      Type  Module   Description
    !     ----------------------------------------------------------------
    !     ----------------------------------------------------------------
    !
    !  6. Error messages :
    !  7. Remarks
    !  8. Structure :
    !  9. Switches :
    !
    !     !/S  Enable subroutine tracing.
    !
    ! 10. Source code :
    !
    !/ ------------------------------------------------------------------- /
#ifdef W3_S
    USE w3servmd, ONLY: strace
#endif
    !
    USE w3odatmd, ONLY: nbi, ndse, isbpi, xbpi, ybpi
    USE w3gdatmd, ONLY: nx, xgrd, ygrd, mapsta, mapfs, mapsf
 
 
    REAL, INTENT(IN)         :: DISTMIN
    LOGICAL, INTENT(INOUT)         :: FLOK
 
    INTEGER                   :: I, J, JMEMO, IS, IX,  N, IX1(NBI)
    REAL                      :: DIST, DIST0
    !
    n = 0
    !
    !1. look for input boundary point index
    ! warning: if land points are included as boundary points to abide by the nest
    ! file, their status should be -2.
    !
    ix1 = 0
    isbpi = 1
    DO ix = 1, nx
      IF (abs(mapsta(1,ix)) .EQ. 2) THEN
        n = n + 1
        IF (n.GT.nbi) THEN
          WRITE(ndse,*) 'Error: boundary node index > NBI ... nest.ww3 file is not consistent with mod_def.ww3'
          stop
        ENDIF
        ix1(n) = ix
#ifdef W3_T
        WRITE(ndse ,*)'ADDING BOUNDARY POINT:',n,ix
#endif
      END IF
    END DO
    !
    !2. Matches the model grid points (where MAPSTA = 2) with the points in nest.ww3
    !   For this, we use the nearest point in the nest file.
    !
    DO i = 1, nbi
      dist0 = huge(1.)
      is = 1
      DO j = 1, n
        dist = (xbpi(i) - xgrd(1,ix1(j)))**2 + (ybpi(i) - ygrd(1,ix1(j)))**2
        IF (dist.LT.dist0) THEN
          is = mapfs(1,ix1(j))
          dist0 = dist
          jmemo = j
        END IF
      END DO
 
      dist0 = sqrt(dist0)
      IF (dist0.LE.distmin) THEN
        isbpi(i) = is
#ifdef W3_T
        WRITE(ndse ,'(A,I6,A,I7,A,I6)') 'MATCHED BOUNDARY POINT:',i,'GRID POINT:', &
             mapsf(is,1),'INDEX IN nest.ww3:', jmemo
#endif
      ELSE
        flok=.true.
      END IF
 
    END DO
 
    IF ( n .NE. nbi) THEN
      WRITE(ndse ,900) n, nbi
      DO j=1,n
        WRITE(6,*) 'THIS POINT HAS MAPSTA=2:',isbpi(j)
      END DO
      isbpi(n+1:nbi)=isbpi(1)
    END IF
 
900 FORMAT (/' *** WAVEWATCH III ERROR IN W3IOBC : '/                &
         '     NUMBER OF MAPSTA=2 DIFFERS FROM NUMBER IN nest.ww3    '/                &
         '     CHECK nest.ww3 AND ww3_grid.inp ',2i8/)

References w3odatmd::isbpi, w3gdatmd::mapfs, w3gdatmd::mapsf, w3gdatmd::mapsta, w3odatmd::nbi, w3odatmd::ndse, w3gdatmd::nx, w3servmd::strace(), w3odatmd::xbpi, w3gdatmd::xgrd, w3odatmd::ybpi, and w3gdatmd::ygrd.

Referenced by w3iobcmd::w3iobc().

Variable Documentation

n_outside_boundary

integer, save w3triamd::n_outside_boundary

Definition at line 114 of file w3triamd.F90.

  INTEGER, SAVE               :: N_OUTSIDE_BOUNDARY

Referenced by readmsh(), and ug_getopenboundary().

outside_boundary

integer, dimension(:), allocatable, save w3triamd::outside_boundary

Definition at line 115 of file w3triamd.F90.

  INTEGER, SAVE, ALLOCATABLE  :: OUTSIDE_BOUNDARY(:)

Referenced by readmsh(), and ug_getopenboundary().