!----- AGPL -------------------------------------------------------------------- ! ! Copyright (C) Stichting Deltares, 2015. ! ! This file is part of Delft3D (D-Flow Flexible Mesh component). ! ! Delft3D is free software: you can redistribute it and/or modify ! it under the terms of the GNU Affero General Public License as ! published by the Free Software Foundation version 3. ! ! Delft3D is distributed in the hope that it will be useful, ! but WITHOUT ANY WARRANTY; without even the implied warranty of ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ! GNU Affero General Public License for more details. ! ! You should have received a copy of the GNU Affero General Public License ! along with Delft3D. If not, see . ! ! contact: delft3d.support@deltares.nl ! Stichting Deltares ! P.O. Box 177 ! 2600 MH Delft, The Netherlands ! ! All indications and logos of, and references to, "Delft3D", ! "D-Flow Flexible Mesh" and "Deltares" are registered trademarks of Stichting ! Deltares, and remain the property of Stichting Deltares. All rights reserved. ! !------------------------------------------------------------------------------- ! $Id: network.f90 42642 2015-10-21 11:34:20Z dam_ar $ ! $HeadURL: https://repos.deltares.nl/repos/ds/trunk/additional/unstruc/src/network.f90 $ module m_netw use network_data use m_alloc implicit none contains subroutine loadNetwork(filename, istat, jadoorladen) use unstruc_netcdf, only : unc_read_net, unc_write_net use unstruc_messages use m_missing implicit none character(*), intent(in) :: filename !< Name of file to be read (in current directory or with full path). integer, intent(out) :: istat !< Return status (0=success) integer, intent(in) :: jadoorladen ! double precision, allocatable, save :: zkold(:) integer :: minp, K, K0, L0, L, NUMKN, NUMLN logical :: jawel character(len=len_trim(filename)) :: netname real :: rmin, rmax CALL CLEARFLOWMODELINPUTS() ! TODO: Zou eigenlijk verplaatst moeten worden [AvD] inquire(file = filename, exist=jawel) if (.not. jawel) then call mess(LEVEL_WARN,'could not open '''//trim(filename)//'''') return end if IF (JADOORLADEN == 0) THEN K0 = 0 L0 = 0 ELSE K0 = NUMK L0 = NUML ENDIF ! New NetCDF net file call unc_read_net(filename, K0, L0, NUMKN, NUMLN, istat) !if (.not. allocated(zkold) ) then ! allocate (zkold(numkn)) ! zkold = zk !else ! do k = 1,numkn ! zk(k) = zk(k) - zkold(k) ! enddo !endif if (istat == 0) then NUMK = K0 + NUMKN NUML = L0 + NUMLN CALL SETNODADM (0) else ! Retry: Original .net files (new and old format) CALL OLDFIL (MINP, filename) CALL REAnet (MINP, istat, jadoorladen) if (istat > 0) then ! reanet yields >0 on success! [AvD] istat = 0 else istat = 1 endif if (istat == 0) then ! Autosave NetCDF form of the network just read. L = index(filename, '.') call unc_write_net(filename(1:L-1)//'_net.nc') call mess(LEVEL_INFO, 'Autosaved NetCDF form of net to ', filename(1:L-1)//'_net.nc') endif call doclose(minp) endif CALL CLOSEWORLD() ! STITCH 0-360 FOR 0-360 GLOBE MODELS netstat = NETSTAT_CELLS_DIRTY end subroutine loadNetwork end module m_netw module m_netstore use network_data implicit none integer, parameter :: maxnodespercell = 6 !< maximum number of nodes per cell integer, save :: maxlinkspernode = 1 !< maximum number of links attached to a node integer, save :: maxnodes = 1 integer, save :: maxlinks = 1 integer, save :: maxcells = 1 type (tnod), allocatable, dimension(:) :: nod_st ! dimension(maxnodes) type (tface), allocatable, dimension(:) :: netcell_st ! dimension(maxcells) double precision, allocatable, dimension(:) :: xk_st, yk_st ! dimension(maxnodes) integer, allocatable, dimension(:) :: nmk_st ! dimension(maxnodes) integer, allocatable, dimension(:) :: nb_st ! dimension(maxnodes) integer, allocatable, dimension(:) :: lnn_st ! dimension(maxlinks) integer, allocatable, dimension(:,:) :: lne_st ! dimension(2,maxlinks) integer, allocatable, dimension(:,:) :: kn_st ! dimension(3,maxlinks) integer :: numnodes, numlinks, numcells ! pointers integer, allocatable, dimension(:) :: ik_st !< nodes, dimension(maxnodes) integer, allocatable, dimension(:) :: iL_st !< links, dimension(maxlinks) integer, allocatable, dimension(:) :: ip_st !< cells, dimension(maxcells) contains !> save network data locally subroutine local_netstore(k) use m_netw use m_alloc implicit none integer, dimension(:), intent(in) :: k !< list of cells integer :: i, j, ii integer :: inode, ilink, icell, M, N do i=1,size(k) N = netcell(k(i))%N if ( N.gt.maxnodespercell ) then call qnerror('local_netstore: array size error', ' ', ' ') return end if end do M=sum(netcell(k)%N) ! determine upper bounds for arrays maxnodes = max(M, maxnodes) maxlinks = max(M, maxlinks) maxcells = max(size(k), maxcells) maxlinkspernode = max(0, maxlinkspernode) do i=1,size(k) icell = k(i) do j=1,netcell(icell)%N inode = netcell(icell)%nod(j) maxlinkspernode = max(nmk(inode), maxlinkspernode) end do end do if ( .not.allocated(ik_st) ) then ! allocation allocate(nod_st(maxnodes), xk_st(maxnodes), yk_st(maxnodes), nmk_st(maxnodes), nb_st(maxnodes), ik_st(maxnodes)) allocate(lnn_st(maxlinks), lne_st(2,maxlinks), kn_st(3,maxlinks), iL_st(maxlinks)) allocate(netcell_st(maxcells), ip_st(maxcells)) do i=1,maxcells allocate(netcell_st(i)%lin(maxnodespercell), netcell_st(i)%nod(maxnodespercell)) end do do i=1,maxnodes allocate(nod_st(i)%lin(maxlinkspernode)) end do else ! memory checks and reallocation if necessary ! maxnodes check if ( size(ik_st).lt.maxnodes ) then deallocate(nod_st) allocate(nod_st(maxnodes)) do i=1,maxnodes allocate(nod_st(i)%lin(maxlinkspernode)) end do call realloc(xk_st, maxnodes) call realloc(yk_st, maxnodes) call realloc(nmk_st, maxnodes) call realloc(nb_st, maxnodes) call realloc(ik_st, maxnodes) end if ! maxlinks check if ( size(iL_st).lt.maxlinks ) then call realloc(lnn_st, maxlinks) call realloc(lne_st, (/2,maxlinks/)) call realloc(kn_st, (/3,maxlinks/)) call realloc(iL_st, maxlinks) end if ! maxcells check if ( size(ip_st).lt.maxcells ) then deallocate(netcell_st) allocate(netcell_st(maxcells)) do i=1,maxcells allocate(netcell_st(i)%lin(maxnodespercell), netcell_st(i)%nod(maxnodespercell)) end do call realloc(ip_st, maxcells) end if ! maxlinkspernode check if ( size(nod_st(1)%lin) .lt. maxlinkspernode ) then do i=1,maxnodes call realloc(nod_st(i)%lin, maxlinkspernode) end do end if end if numnodes = 0 numlinks = 0 numcells = size(k) ip_st(1:numcells) = k do i=1,numcells icell = ip_st(i) N = netcell(icell)%N if ( icell.eq.5005 ) then continue end if ! store cell data netcell_st(i)%lin(1:N) = netcell(icell)%lin(1:N) netcell_st(i)%nod(1:N) = netcell(icell)%nod(1:N) netcell_st(i)%N = N ! store node and link data ! use nmk<-100 to mask nodes ! use lnn<-100 to mask links do j=1,netcell(icell)%N inode = netcell(icell)%nod(j) ilink = netcell(icell)%lin(j) if ( nmk(inode).gt.maxlinkspernode ) then call qnerror('local_netstore: array size error', ' ', ' ') return end if if ( nmk(inode).ge.0 ) then numnodes = numnodes+1 ik_st(numnodes) = inode nod_st(numnodes)%lin(1:nmk(inode)) = nod(inode)%lin(1:nmk(inode)) nmk_st(numnodes) = nmk(inode) nb_st(numnodes) = nb(inode) xk_st(numnodes) = xk(inode) yk_st(numnodes) = yk(inode) nmk(inode) = -nmk(inode) - 100 end if if ( lnn(ilink).ge.0 ) then numlinks = numlinks+1 iL_st(numlinks) = ilink lnn_st(numlinks) = lnn(ilink) lne_st(:,numlinks) = lne(:,ilink) kn_st(:,numlinks) = kn(:,ilink) lnn(ilink) = -lnn(ilink) - 100 end if end do end do ! reset mask do i=1,numnodes nmk(ik_st(i)) = -nmk(ik_st(i)) - 100 end do do i=1,numlinks lnn(iL_st(i)) = -lnn(iL_st(i)) - 100 end do return end subroutine local_netstore !> restore network data locally subroutine local_netrestore() use m_netw use m_alloc implicit none integer :: i, N integer :: icell, iLink, inode do i=1,numcells icell = ip_st(i) N = netcell_st(i)%N netcell(icell)%N = N call realloc(netcell(icell)%nod, N) call realloc(netcell(icell)%lin, N) netcell(icell)%nod(1:N) = netcell_st(i)%nod(1:N) netcell(icell)%lin(1:N) = netcell_st(i)%lin(1:N) end do do i=1,numlinks iLink = iL_st(i) lnn(ilink) = lnn_st(i) lne(:,ilink) = lne_st(:,i) kn(:,ilink) = kn_st(:,i) end do do i=1,numnodes inode = ik_st(i) N = nmk_st(i) call realloc(nod(inode)%lin, N) nod(inode)%lin(1:N) = nod_st(i)%lin(1:N) nmk(inode) = nmk_st(i) nb(inode) = nb_st(i) xk(inode) = xk_st(i) yk(inode) = yk_st(i) end do return end subroutine local_netrestore end module m_netstore