anchor.diff:  c
         integer icflip(4)
         data icflip/1,3,2,4/
  -      integer kwb,kwc,kwfb(maxflts),kwfc(maxflts)
  +      integer kwb,kwc,kww,kwfb(maxflts),kwfc(maxflts)
         data kwb/12/
         data kwc/13/
  +      data kww/12/
         data kwfb/12,13,14,15,16,17/
         data kwfc/18,19,20,21,22,23/
         integer izero
         data izero/0/
  +      double precision zero
  +      data zero/0.0d0/
   c
   c...  parameters read from parameter file
   c
  -      integer nbc,iconopt
  -      integer ibcode(maxbnds),ibc(nsd,maxbnds),isn(nsd)
  +      integer nbc,iconopt,nwsets
  +      integer ibcode(maxbnds),ibc(nsd,maxbnds),isn(nsd),iwcode(maxwsets)
  +      integer iwdir(maxwsets),modew(maxwsets)
         integer iftype(maxflts),nmatf(2,maxflts),ifmat(2,maxflts,maxmatf)
         double precision bc(nsd,maxbnds),fsplit(nsd,2,maxflts)
  -      double precision cscale
  +      double precision rhow(maxwsets),gw(maxwsets),cscale
   c
   c...  filenames
   c
         character fileroot*200,pfile*200,ifile*200,nfile*200,cfile*200
  -      character bcfile*200,bccfile*200
  +      character bcfile*200,bccfile*200,wbcfile*200
         character fbcfile*200,fbccfile*200
   c
   c...  parameters and variables read from netgen file
  
positive.code: c
  c...  external routines
  c
        double precision tetcmp
        integer nnblnk,nchar
        external nnblnk,nchar,tetcmp
  c
  c...  local variables
  c
        double precision det,sgn,xl(nsd,nen),xtmp(nsd)
        integer iadjf(maxflts*maxfnodes,maxnbrs),inodef(maxflts*maxfnodes)
        integer ientmp(nen),nelsf(maxflts*maxfnodes),ifltnode(maxnodes)
        integer ifhist(maxflts),idir(nsd),ibcnode(maxnodes,maxbnds)
        integer nbcnodes(maxbnds)
        integer i,j,k,l,jj,i1,i2,j1,j2,nenl,nsdl,nfltnodes,kk
        integer iel,nflip
        character cstring*15,cbound(maxbnds)*3,cfnum*2
        data cstring/"coord_units = m"/
        data cbound/
       & "b01","b02","b03","b04","b05","b06","b07","b08","b09","b10",
       & "b11","b12","b13","b14","b15","b16","b17","b18","b19","b20",
       & "b21","b22","b23","b24","b25","b26","b27","b28","b29","b30"/
  c
        nenl=nen
        nsdl=nsd
  
negative.code:             end if
            end do
          end if
          do j=1,nwsets
            if(ibcfac.eq.iwcode(j)) then
              area=acomp(x,ibcvert,nsd,numnp,iwdir(j))
              do k=1,3
                iwbcnode(ibcvert(k),j)=ibcfac
                wink(ibcvert(k),j)=wink(ibcvert(k),j)+
       &         area*rhow(j)*gw(j)/3.0d0
              end do
            end if
          end do
        end do
        close(kr)
        nfltnodes=0
        do i=1,numnp
          if(ifltnode(i).eq.1) then
  