vtf-logo

fsi/sfc-amroc/WaterBlastFracture/src/combl.f

c     =====================================================
      subroutine combl()
c     =====================================================
c     
c     Create and initialize application specific common-blocks.
c
c     Copyright (C) 2003-2007 California Institute of Technology
c     Ralf Deiterding, ralf@amroc.net
c
      implicit double precision (a-h,o-z)
      include  "cuser.i"
c
      parameter( lin = 5, lmechout = 11, lrin = 13 )
      character *16 cwork
      dimension cwork(3)
      character *256 filename
c
      RU = 0.83140000d+01
      PA = 0.10132500d+06
c
      cwork(1)= 'Water'
      cwork(2)= 'Air'
      Wk(1) = 18.015d-3
      Wk(2) = 29.08d-3
      g(1) = 7.415d0
      g(2) = 1.4d0
      pinf(1) = 296.2d6
      pinf(2) = 0.d0
c
      cwork(3)= 'Pinf'
      Wkhelp = 1.d0
c
      open(unit=lmechout, status='unknown', form='formatted', 
     &     file='chem.dat')
      write (lmechout,400) RU
      write (lmechout,401) PA
      write (lmechout,402) (k,cwork(k),k=1,3)        
      write (lmechout,403) (k,Wk(k),k=1,2)  
      write (lmechout,403) 3,Wkhelp       
      close (lmechout)
c
 400  format('RU       ',e16.8)
 401  format('PA       ',e16.8)
 402  format('Sp(',i2.2,')     ',a16)
 403  format('W(',i2.2,')    ',e16.8)
c
      open(unit=lin,status='old',form='formatted',file='init.dat')
      read (lin, *) rd, al, Np
      read (lin, *) v0, cm, Ninit
      read (lin, *) Nread, filename, sshift
      close (lin)
c
      cw = 1482.d0
      rhow = 1.d3
c 
      rhoo = 1.d0
      uo   = 0.d0
      po   = PA
      as   = 0.d0
c
      rhoi = rhow
      ui   = 0.d0
      pi   = PA
c
      xm = al
      vm = -v0
      area = 4.d0*datan(1.d0)*rd**2
c
c     # Use p0 and theta
      if (Ninit.eq.1) then
         vm = -v0/(rhow*cw)
         cm = cm*area*rhow*cw
      endif
c
      if (Nread.gt.0) then
         do k=1,256
            if (filename(k:k).eq.' ') goto 100
         enddo
 100     Nread = Nread + 2 - 1
         open(unit=lrin, status='old', form='formatted', 
     &        file=filename(1:k-1))
         read (lrin,*) xm, vm
         do i=2,Nread-1
            read (lrin,*) xtab(i),(qtab(k,i),k=1,3)
            xtab(i) = xtab(i)-sshift
         enddo
         close (lrin)
         xtab(1) = xtab(2) - 0.5d0*(xtab(3)-xtab(2))
         alpha = (xtab(1)-xtab(2))/(xtab(3)-xtab(2))
         do k=1,3
            qtab(k,1) = (1.d0-alpha)*qtab(k,2) + alpha*qtab(k,3)
         enddo
c
         xtab(Nread) = xtab(Nread-1)+0.5d0*(xtab(Nread-1)-xtab(Nread-2))
         alpha = (xtab(Nread)-xtab(Nread-1))/
     &        (xtab(Nread-2)-xtab(Nread-1))
         do k=1,3
            qtab(k,Nread) = (1.d0-alpha)*qtab(k,Nread-1) + 
     &           alpha*qtab(k,Nread-2)
         enddo
      else
         write (6,*) 'Boundary velocity = ',vm
         write (6,*) 'Mass = ',cm
         write (6,*) 'Mass per unit area = ',cm/area
         write (6,*) 'Theta = ',cm/(area*rhow*cw)
         write (6,*) 'p0 = ',-vm*rhow*cw
      endif
c
      return
      end

c     =====================================================
      subroutine cblread(xmt,vmt,cmt,rdt,npt)
c     =====================================================
c     
      implicit double precision (a-h,o-z)
c
      include  "cuser.i"
c
      xmt = xm
      vmt = vm
      cmt = cm
      rdt = rd
      npt = Np
c
      return
      end

c     =====================================================
      subroutine cblwrite(xmt,vmt,cmt,rdt,npt)
c     =====================================================
c     
      implicit double precision (a-h,o-z)
c
      include  "cuser.i"
c
      xm = xmt
      vm = vmt
      cm = cmt 
      rd = rdt
      Np = npt
c
      return
      end